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ABSTRACT 

The hot Jupiter HD 189733 b is the most extensively observed exoplanet. Its atmo- 
sphere has been detected and characterised in transmission and eclipse spectroscopy, 
and its phase curve measured at several wavelengths. This paper brings together the 
results of our campaign to obtain the complete transmission spectrum of the atmo- 
sphere of this planet from UV to infrared with the Hubble Space Telescope, using the 
STIS, ACS and WFC3 instruments. We provide a new tabulation of the transmis- 
sion spectrum across the entire visible and infrared range. The radius ratio in each 
wavelength band was re-derived, where necessary, to ensure a consistent treatment 
of the bulk transit parameters and stellar limb-darkening. Special care was taken to 
correct for, and derive realistic estimates of the uncertainties due to, both occulted 
and unocculted star spots. 

The combined spectrum is very different from the predictions of cloud-free models 
for hot Jupiters: it is dominated by Rayleigh scattering over the whole visible and near 
infrared range, the only detected features being narrow sodium and potassium lines. 
We interpret this as the signature of a haze of condensate grains extending over at 
least five scale heights. We show that a dust-dominated atmosphere could also explain 
several puzzling features of the emission spectrum and phase curves, including the 
large amplitude of the phase curve at 3.6 |xm, the small hot-spot longitude shift and 
the hot mid-infrared emission spectrum. We discuss possible compositions and derive 
some first-order estimates for the properties of the putative condensate haze/clouds. 
We finish by speculating that the dichotomy between the two observationally defined 
classes of hot Jupiter atmospheres, of which HD 189733 b and HD 209458 b are the 
prototypes, might not be whether they possess a temperature inversion, but whether 
they are clear or dusty. We also consider the possibility of a continuum of cloud 
properties between hot Jupiters, young Jupiters and L-type brown dwarfs. 

Key words: planetary systems - stars: individual (HD189733) - planets and satel- 
lites: atmospheres - techniques: spectroscopic 



1 INTRODUCTION 

Transiting exoplanets are giving us our first glimpse into 
planetary atmospheres beyond the Solar System. Observa- 
tional biases make short-orbit gas giants (hot Jupiters) most 
accessible to observations (see [Haswell 2010 , for a textbook- 
level introduction). The planet HD 189733b occupies a place 
of choice among the observable targets, because of its bright 
host star and large planet-to-star radius ratio. HD 189733b 

* E-mail: fpont@astro.ex.ac.uk 



has been observed abundantly with both the Hubble and 
Spitzer space telescopes in order to characterise its atmo- 
sphere. It epitomises both the promises and challenges of 
these observations. 

1.1 The hot Jupiter HD 189733b 

HD 189733b is a 1.15 Mjup gas giant, orbiting its K-dwarf 
host star in 2.2 days. Due to the relatively small size of the 
star, a 0.77 Rq K dwarf, the transit signal is deep (2.4%). 
Three types of observations can be used to characterise 
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the atmosphere of the planet: the transmission spectrum 
along the limb during transits, the day-side emission spec- 
trum during secondary eclipses, and the day-night tempera- 
ture contrast with the orbital phase curve. One atmospheric 
scale height projects to about 1x10"* in transmission spec- 
troscopy. The secondary eclipse is 3xl0~^ deep at 8 vim, 
and the phase curve amplitude is 1 x 10~^ at 4 fim. Data 
on the atmosphere of HD 189733b using these three meth- 
ods have been gathered using the ACS, STS, NICMOS and 
WFC3 instruments on the Hubble Space Telescope ("HST") 
and the IRAC and MIPS instruments on the Spitzer Space 
Telescope ("^piteeT^'Q 



1.2 Expectations from atmosphere models 

The irradiation temperature (equilibrium temperature at 
the sub-stellar point) of HD 189733b is 1650 K, making 
the zero-albedo equilibrium temperature 1170 to 1390 K 
depending on the efficiency of heat transfer from the day 
side to the night side. At these temperatures, hot jupiter 
atmosphere models predict a clear, dark atmosphere, dom- 
inated by absorption by neutral sodium and potassium in 
the visible, and molecular bands of water, ammonia, car- 
bon monoxide, carbon dioxide and methane in the infrared 



(Fortney et al. 20081. In the blue and near UV, scattering 
by hydrogen molecules is expected to become the dominant 
source of opacity. The temperature is too cold to sustain ti- 
tanium oxide and vanadium oxide vapours, that can provide 
high opacities in the visible for hotter atmospheres. 

The dominant source of heat in the atmosphere is the 
stellar irradiation, which far exceeds the internal heat leak- 
age from the interior of the planet. As a result of the in- 
jection of heat from above, the atmosphere is expected to 
be vertically stable, which radiative transfer being the dom- 
inant mechanism of vertical energy exchange. For this rea- 
son, although the atmospheric temperature of the planet is 
close to the condensation temperature of several abundant 
components, including silicates and iron, atmosphere struc- 
ture models predict an atmosphere free of condensates as 
a result of vertical stability. The possibility of clouds was 
considered from the earliest models though ( [Seager &: Sas-| 
selovj 2000 ) , for instance due to the mixing associated with 



the vigorous day-night atmospheric circulation. According 
to models, condensates would make the spectral features 
weaker, or mask some of them, depending on the height of 



the cloud deck ( 


Marley et al.|1999,; Seager & Sasselov|2000 


|Sudarsky et al.| 


20031). Smaller grains or partly transparent 



hazes can also imprint some features on the spectrum, such 
as Rayleigh scattering in the blue or silicate absorption fea- 
tures in the mid-infrared. 

Photochemistry may also play a role in producing con- 
densates, for instance sulphur or carbon compounds, that 
could affect the spectral signatures of the planetary atmo- 
sphere ( |Zahnle et al.|2009b I. 



1.3 Observation of the atmosphere of HD 189733b 

The transit spectrum of HD 189733b has been measured 
from 300 nm to 1 \i.m with the STIS and ACS instruments 
aboard HST, from 1 to 3 jxm with HST's WFC3 and NIC- 
MOS, and in five infrared passbands between 3.6 and 24 \im 
with IRAC and MIPS on Spttzer. The emission spectrum 
has been measured in the five IRAC and MIPS passbands 



and from 8 to 13 ^m with Spitzer's IRS (Agol et al. 2010 



Grillmair et al. 2008 Desert et al. 2009 20111. The phase 



curve along the orbit has been measured in the 3.6, 4.5 and 8 



[i.m channels of IRAC, and around 24 ii.m with MIPS ( Knut- 
son et al.|2007a 2012 \. The signature of Rayleigh scattering 



i Pont et al.|2008[ |Sing et al.||2011|) and the sodium doublet 



(Redfield et al. 2008 Huitson et al. 20121 have been de- 



tected in the UV and visible. Claims of detection of molec- 
ular features in the infrared ( [Swain et al.|[2008l ) have not 
been confirmed by subsequent measurements and analyses 



fSing et al.|2009| [Desert et al.|2009l|2011|[Glbson et al.|20"TT 
|2012bJ 

Inferring atmospheric properties from the present ob- 
servations is a challenge. The signals remain small even by 
the standard of space observations, and the uncertainties 
are generally dominated by complex instrumental effects. 
It has not been uncommon for new observations, and even 
new analysis of existing observations, to contradict earlier 
results. In the infrared, the wavelength coverage is sparse 
(imposed by the Spttzer passbands in space and the water 
absorption gaps from the ground). With a handful of pass- 
bands for as many molecules, interpretation often admits 
several possibilities. 

In the emission spectrum, the planetary fiux data in 
Spitzer passbands were compatible with a temperature pro- 
file decreasing with height and the presence of water and 
other molecules. However, new phase curves at 3.6 and 
4.5 [im and a re-analysis of the previous Spitzer data by 



Knutson et al. (20121 have also overhauled previous results 
on the emission spectrum. The transit depths from |Char-| 
[bonneau et al.| ( [2008| ), crucial to constrain the temperature 
profile and water abundance, turn out to be incorrect by 
more than 3 sigma, due to the complex Spitzer instrumen- 
tal systematics (the "ramp effect" and "pixel phase effect" 
affecting high-accuracy fiux measurements). The amplitude 
of the phase curves and depth of the eclipses at 3.6 and 
4.5 \im are no longer amenable to a simple interpretation. 



Knutson et al. (20121 mention non-equilibrium chemistry 
and slowed day-night recirculation as possible explanations 
of the anomalies in the data. 

Overall, the current interpretation of the ensemble data 
is of an atmosphere corresponding to dust-free models with 
molecular absorption and without a stratospheric temper- 
ature inversion. A thin layer of haze at high altitude pro- 
duces the Rayleigh scattering signature in the visible in the 
transmission spectrum, but because the grains are small, 
it does not affect the emission spectrum, the transmission 
spectrum in the infrared or the global energy budget of the 
atmosphere. 



^ Advanced Camera for Surveys, Near Infrared Camera and 
Multi Object Spectrometer, Space Telescope Imaging Spectro- 
graph, Wide Field Camera 3, InfraRed Array Camera, Multiband 
Imaging Photometer for Spitzer 



1.4 HST programme GO-11740 and present paper 

Our HST programme with STIS and WFC3 (GO-11740) 
was designed to put this interpretation to the test. This pro- 
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gramme was aimed at filling the wavelength gaps in trans- 
mission spectroscopy, providing a continuous coverage from 
the UV to near infrared to address the discrepancies be- 
tween the different datasets and lift some ambiguities in 
interpretations. It provides the observational basis of the 
present paper. The explanation in terms of a thin layer of 
haze at high altitude requires the transmission spectrum to 
evolve sharply from a featureless slope to deep molecular 
features in the 1-2 (i.m interval between the ACS and NIC- 
MOS measurements. The STIS and WFC3 measurements 
have shown the featureless slope to extend beyond the range 
probed initially by ACS, suggesting that Rayleigh scattering 
dominates the transmission spectrum over at least five atmo- 
spheric scale heights. The finding that Rayleigh scattering 
dominates the transmission spectrum continuously from the 
UV to the near infrared opens the possibility that the opti- 
cally thick haze layer is more than a marginal phenomenon 
affecting only the very low-pressure parts of the atmosphere, 
and may affect the radiative transfer and spectrum through- 
out the atmosphere. 

In this paper, we present the global results of the HST 
GO-11740 programme to obtain a complete transmission 
spectrum from UV to near-IR with HST. Results for indi- 



vidual instruments were presented in Gibson et al. (2012a 
and |Sing et al.| ( |2011[ ). We concentrate here on combining the 
different wavelength ranges into a single coherent scale. This 
is rendered difficult primarily because HD 189733 presents 
variability at the 1-2% level due to star spots, so that mea- 
suring the transmission spectrum of the planet implies first 
disentangling the effect of the star spots, which has a similar 
amplitude. We then consider the whole ensemble of data on 
the atmosphere of HD 189733b, including the key update 
by Knutson et al. ( 2012[ ) on secondary eclipses and phase 
curves. We explore the possibility that the dust detected in 
the transmission spectrum has a large impact on the atmo- 
sphere as a whole. 

Section [2] presents the collection and modelling of the 
transit measurements of the GO-11740 programme and 
other programmes pertaining to the transmission spectrum. 
Section |3] tackles the crucial but somewhat technical issue 
of accounting for star spots on the parent star. Section [1] 
presents the global UV-to-infrared transmission spectrum 
obtained by combining the datasets in a consistent manner. 
Section [5] discusses the influence of haze and clouds in the 
atmosphere data as a whole, including emission and phase 
curve information. In Section [6j we explore the possibility 
that of the atmosphere of HD 189733b is dominated by con- 
densates, and examine the observable consequences. 



2 ATMOSPHERE TRANSMISSION 
SPECTRUM: COLLECTING AND 
MODELLING TRANSIT OBSERVATIONS 

The transmission spectrum of the atmosphere of a transiting 
planet can be measured by monitoring the transit in spec- 
troscopy. In principle the method is simple: since the depth 
of the transit is proportional to the square of the planet- 
to-star radius ratio, measuring the depth of the transit at 
different wavelengths records the opacity of the atmosphere 
to grazing star light at those wavelengths, yielding a trans- 
mission spectrum. 



As outlined in the Introduction, the transit depth of the 
HD 189733 system has been measured precisely over many 
wavelengths from 0.3 to 24 ii.m. On paper, the potential of 
these measurements is impressive. The expected amplitude 
of the atmospheric features is about 5 • 10"* as a fraction of 
the stellar flux, and the photon-shot noise integrated over 
the whole transit is smaller than 1 ■ 10"'* at all wavelengths 
up to 8 \im. Such measurements are able to constrain the 
atmospheric structure and composition in significant detail. 

However, in practice, transit spectroscopy of 
HD 189733b has had to contend with two obstacles. 
The first is the ubiquity of complex instrumental effects 
at and above the 10~* level in flux measurements. With 
infrared detectors in particular, these instrumental sys- 
tematics have generally turned out to be impossible to 
correct perfectly, so that unknown uncertainties similar 
in amplitude to the expected signals remain associated to 
most measurements. The second obstacle is the fact that 
HD 189733 is an active, spotted star, which introduces an 
additional signal in the relation between transit depth and 
planet size. 



These issues are discussed in Pont et al. (20081, Desert 



et al. (20091, Sing et al. (20111 and Gibson et al. (20111 



They prevent an easy comparison of the transit radius in 
different wavelength ranges, so that producing a global es- 
timate of the UV-to-IR transmission spectrum requires a 
detailed combined analysis of all data sets. This analysis 
makes up the rest of the present Section. 

In contemplating the possibility of measuring a global 
transmission spectrum for HD 189733b from measurements 
widely separated in time, we make the assumption that in- 
trinsic variations in the atmosphere of the planet over time 
are smaller than the present observational uncertainties. 
This assumption is discussed in Paragraph |3.4| 



2.1 Transmission spectroscopy datasets 

Table |2] gives the journal of the space-based observations of 
the transit of HD 189733b considered in this paper. The 
primary focus of this study is addressed by our spectro- 
scopic measurements with the STIS, ACS and WFC3 instru- 
ments on HST, to provide a continuous transmission spec- 
trum from the near UV (300 nm) to the near IR (2 \im). 
This covers the interval between high-altitude atomic lines 
in the UV and molecular lines near the photosphere in the 
infrared, encompassing the part of the atmosphere involved 
in the bulk of the energy exchange between the stellar irra- 
diation and the planetary atmosphere. We also include spec- 
troscopic and photometric measurements with the NICMOS 
camera on HST, as well as broadband photometric measure- 
ments with the IRAC and MIPS instruments on Spitzer. 

The individual datasets, presented separately in the 
publications mentioned in the previous section, are briefly 
summarised below. 



2.1.1 STIS (HST) 

Five transits of HD 189733 were observed with STIS on the 
refurbished HST in two contexts, at low resolution from the 
present programme, and at medium resolution around the 
sodium doublet (GO-11111). The results of the flrst obser- 



vations are published in Sing et al. (2011 1 and the second in 
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Huitson et al. (20121. The broad-band data showed a fea- 
tureless spectrum over the 300-600 nm range, with absorp- 
tion rising bluewards. The medium-resolution data resolved 
both components of the sodium doublet, showing strong nar- 
row cores and an absence of broad wings for the sodium 
feature. 



2.1.2 ACS (HST) 

Three transits were observed with the ACS camera on HST, 
using Grism G800L (0.6-1 \ira) on Ma y 22, May 26 and J uly 
14, 2006. The result s were published in|Pont et al] ( |2007[ ) for 
the time series and Pont et al. ( 2008 1 for the transmission 
spectrum. 



One important issue for the present study is the mea- 
surement of potassium absorption. A very large line with 
extended wings is expected in a clear atmosphere from the- 
oretical models, but was not seen in the data. The ACS 
measurements do not have the resolution to detect the nar- 
row core of the potassium line, and is only sensitive to its 
total equivalent width over a larger wavelength band. 

We find excess absorption located close to the expected 
position of the core of the potassium doublet, but the reso- 
lution is too low and the noise on pixel-by-pixel decorrela- 
tion of systematics too high to give a measurement of the 
width of the detection. We therefore fit the amplitude of a 
potassium line model to the pixel-by-pixel data, with the 



model used in Huitson et al. ( 2012 1 to fit the potassium line 



and surrounding Rayleigh scattering slope, which is based 
upon the analytic expression from [Lecavelier Pes Etangs| 
et al.'f2008). As with the sodium feature, for the potassium 



cross-section we neglected pressure broadening due to the 
observed lack of broad line wings. We fixed the abundance 
to an arbitrary level, which is unconstrained as the haze has 
an unknown composition, and fit for the baseline altitude of 
a 900 A spectral region surrounding the sodium doublet, as 
well as the model temperature, assumed to be isothermal. A 
900 A width was chosen such that the baseline altitude was 
well measured, while also ensuring that the best-fit model 
temperature was not sensitive to the broadband Rayleigh 
slope. The result is shown on Fig. [l] The potassium feature 
in the ACS data is 2.5-(j significant. The model temperature 
scales to fit the amplitude of the potassium feature, with a 
best-fit value of 1800±720 K. 

Thus, while the lack of high-resolution data around the 
potassium line prevents very specific constraints, the best- 
fit solution is compatible with the presence of the line to 
the same level as the sodium line or higher, i.e. with a sig- 
nificant core but no broad wings. The possible presence of 
potassium in the ACS data is in contrast to the findings 



of Jensen et al. (20111 who did not detect potassium us- 
ing high-resolution ground-based spectroscopy. However, the 
different resolutions, wavelength range, and planetary alti- 
tudes probed combined with the low signal-to-noise of both 
studies makes a direct comparison difficult. 



2.1.3 WFC3 (HST) 

Two transits were observed with the WFC3 camera on the 
refurbished HST as part of our programme GO-11740, cov- 
ering the wavelength interval 1.1-1.7 \im. The results are 




7000 7200 7400 7600 7800 8000 8200 8400 
Wovelength (A) 



Figure 1. ACS transmission spectrum of HD 189733b in a 900 A 
region surrounding the potassium doublet. Also plotted is a best- 
fit model (both unbinned and binned to the data) which includes 
Rayleigh scattering and sodium absorption. 



presented in Gibson et al. (2012a I. Due to saturation and 



non-linearity affecting the brightest (central) pixels of the 
spectrum, light curves were extracted from the blue and red 
ends of the spectra only, corresponding to wavelength ranges 
of 1.099-1.168 and 1.521-1.693 ^.m for the first visit, and 
1.082-1.128 and 1.514-1.671 ^m for the second. To account 
for instrumental systematics, the light curves were fitted us- 
ing a Gaussian-process model whilst simultaneously fitting 
for the transit parameters. 



2.1.4 NICMOS (HST) 

Transit spectroscopy data for HD189733b has been gathered 
by NICMOS in two modes: spectroscopic and photometric. 
The first is presented in Swain et al. (20081, the second in 



Sing et al. j2009 ) . The two datasets reach incompatible con 



and Gibson et al. (2012b I 



reanalyse the spectroscopic data, and discuss this issue. The 
reader is referred to these papers for details. We adopt the 



method of Gibson et al. 



(2012b I to calculate the uncertain- 
ties of the spectroscopic data. These use a Gaussian pro- 
cess to model the instrumental systematics, which avoids 
the restrictive assumption of linear basis models and per- 
forms the inference in a Bayesian way, therefore mitigating 
against over-interpreting systematics. 



2.1.5 IRAC (Spitzer) 



We use Desert et al. ( 2009 1, Desert et al. ( 2011 1 and Knutson 



et al. (20121 for measurements in the 3.6-[i.m, 4.5-^J,m and 



5.8 n-m channels of IRAC and with IRAS on Spitzer, and 



Agol et al. (20101 for the 8-micron IRAC channel. 



2.2 Self-consistent transit modelling 

Three issues must be addressed in order to place all the 
observations of the radius of HD 189733b on a common scale 
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Transit central time 
Orbital period 
Impact parameter 
Tangential velocity 



2454279.436714 [BJD] 
2.21857567 days 

0.6631 
25.125 R* days-i 



Table 1. Orbital parameters adopted for the re-analysis in this 
paper, from'Agol et al.l(j2010|: timing of transit, orbital period, 
impact parameter, tangential velocity during the transit. 



Instrument 


Setting 


Number of 


Wavelength 






transits used 


range 


[nm] 


STIS 


G430L 


2 


290 - 


570 


STIS 


G750M 


3 


581 - 


638 


ACS 


G800L 


1 


550 - 


1050 


NICMOS 


G206 


1 


1200 - 


2600 


NICMOS 


F166N 


2 


1649 - 


1667 


NICMOS 


F187N 


2 


1864 - 


1874 


WFC3 


G141 


2 


1083 - 


1693 


IRAC 


IDll 


3 


3200 - 


4000 


IRAC 


2D09 


2 


4000 - 


5000 


IRAC 


3D09 


1 


5100 - 


6500 


IRAC 


4D09 


7 


6400 - 


9300 


MIPS 




1 


19500 - 


28500 



Table 2. List of the transit observations of HD 189733 from space 
used in the present study. 



2.2.2 Limb- darkening 

Stellar limb darkening modifies the shape of a planetary 
transit, and therefore afTects the radius determination. This 
is important when measuring the transmission spectrum, be- 
cause limb darkening has a strong wavelength dependence. 
The effect of limb darkening on the measurement of the 
transmission spectrum gets stronger towards shorter wave- 
lengths. It is negligible compared to other sources of uncer- 
tainty in the infrared, but becomes significant in the visible 
and UV. Fortunately, the HST data on HD 189733b is ac- 
curate enough to provide a constraint on limb darkening 
from the data itself, independently of stellar models. Some 
slight discrepancies between theoretical and observed limb 
darkening have been observed in STIS data for HD 209458b 
( Knutson et al.||2007b[ ). [Hayek et al] ( |2012[ ) have computed 
improved limb darkening coefficients from 3-D stellar atmo- 
sphere models, and find that these can capture the behaviour 
of the STIS data for HD 209458b b etter than the 1- D stellar 
atmosphere models used in ,Knutson et al.| ( |2007b] >. For the 
temperature of HD 189733b, the difference between the pre- 
diction of 1-D and 3D models is lower, and both can account 
for the observed shapes of the transits with ACS and STIS. 
Here we use the new set of limb darkening coefficients from 



Hayek et al. (20121 



The remaining uncertainties on limb darkening coeffi- 
cients is not a significant source of error for our combined 
transmission spectra. 



and build a transmission spectrum over a wide wavelength 
range: the orbital elements, the stellar limb darkening, and 
the effect of stellar variability and star spots. 



2.2.1 Orbital parameters 

These measurements of radius ratio used to construct the 
transmission spectrum are determined as part of a model 
transit light curve fit, that also includes other parameters of 
the system, namely the orbital ephemerides, the stellar den- 
sity and the orbital inclination. Any offset in these quantities 
will result in a difference in the inferred radius ratio. This 
difference will be slightly wavelength-dependent, first, be- 
cause of the corresponding changes in the effect of stellar 
limb darkening, and second, because different parts of the 
transit were measured with different instruments, so that 
the covariance between the transit depth and other orbital 
parameters will vary. Because different values for the system 
parameters have been used in separate studies, the resulting 
planet radius values cannot be compared directly. 

Fortunately, in the case of HD 189733b, the orbital pa- 
rameters are known with extreme precision. This is there- 
fore not a significant source of uncertainty in the final re- 
sults. The best measurements of the orbital parameters come 
from the monitoring of 14 transits and eclipses at 8-|J.m with 
Spitzer ( |Agol et al.|2010[ see Table [TJ. We re-calculated the 
radius ratio for all data sets using these ephemerides. When 
a linear approximation was sufficient, we corrected the val- 
ues calculated with other parameters using the first-order 
expressions of the relation between orbital parameters and 
transit shape given in Carter et al. |2008). Otherwise we 



re-fitted the original data with the reference system param- 
eters. 



2.2.3 Instrumental effects 

In most of our data sets, the amplitude of residual instru- 
mental systematics is similar to the signal expected from the 
planetary atmosphere. A careful consideration of the correc- 
tion of instrument effects is therefore essential. The discus- 
sion of these corrections occupy a large part of the literature 
cited above, and we will only make a brief overview here. 

Shortwards of 1 \im (i.e. STIS and ACS) the instrumen- 
tal systematics are relatively well-behaved. They mainly re- 
sult from the pointing drift and thermal expansion of the 
telescope, causing slight movements of the spectrum on the 
detector and changes of focus. The telescope pointing and 
focus can be measured precisely on the data itself by mon- 
itoring the position and width of the spectrum on the de- 
tector. The correlation is easy to trace because of the very 
high signal-to-noise ratio of the fiux measurements. As a 
result, in the case HD 189733, the uncertainties on these 
instrumental corrections are smaller than those due to the 
activity of the host star (see Section [3|, and than the signal 
of the atmosphere of the planet. 

Further in the infrared, instrumental systematics be- 
come the dominant source of uncertainty. The NICMOS and 
WFC3 data show fluctuations in the flux measurements that 
are thought to be correlated to changes in the temperature of 
the detector, and this temperature cannot be tracked as pre- 
cisely as the pointing and focus changes. The Spitzer mea- 
surements are affected by a ramp effect and a pixel phase 
effect, that have been abundantly discussed in the relevant 
literature and limit the maximum accuracy of transit and 
eclipse measurements (none of these instruments had been 
originally designed for high-accuracy time series photometry 
of bright sources). 
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Among the infrared datasets for HD 189733b, the NIC- 
MOS photometric measurements are the ones least beset by 
these problems. Since they are passband-integrated photo- 
metric measurements rather than spectra, the instrumen- 
tal effects have no wavelength-dependent component. The 
drawback is that the observations in different wavelengths 
are not simultaneous, rendering transmission spectrum mea- 
surements degenerate with changes in the star or in the 
planetary atmosphere. The NICMOS grism measurements 
are affected by time- and wavelength-dependent systemat- 
ics, that are difficult to correct entirely because their rela- 
tion with measured instrumental parameters is much looser 
than for ASC and STIS. Gaussian-process models of these 
effects show that linear decorrelation is not sufficient, and 
that non-linear relations between the measured instrumen- 
tal parameters and their effect on the measured fluxes must 
be considered (see Gibson et al.|20rT 2012b 1. 

In Spitzer measurements, the instrumental effects are 
generally larger than the atmospheric signatures. This has 
been a severe limitation in the context of the measurement of 
secondary eclipses and day-side surface brightness for many 
exoplanets. A much better understanding of these effects 
has now been reached (see Desert et al.|2009 2011 Knutson 
et al.|2012[ for HD 189733). There are two effects, a gradual 
increase of the detector's sensitivity during the exposure, 
strongest in the shorter-wavelength channels, and a depen- 
dence on pointing, causing flux variations when the satellite 
pointing varies, because of intra-pixel sensitivity variations 
in the detector. 



2.2.4 Flux from the planet 

When converting transit depths into planet-to-star radius 
ratio, we consider the contribution of the emission from the 
night side of the planet to the total flux. We use a T = 900 K 
blackbody distribution for the planetary emission (appropri- 
ate for the night side of the planet seen during transit), and a 
T = 5000 K blackbody for the stellar emission. The relation 
between the observed and true radius ratio r is: 



2 ^T l 



B 



= n..uc(i + r'-^) 

This correction amounts to 1.2- 10~* on the radius ratio 
at 8 \xm, 0.27 • 10~* at 3.6 ^im. It is negligible in the visi- 
ble, and smaller than the other sources of uncertainty at all 
wavelengths. 



3 IMPACT OF STELLAR VARIABILITY AND 
STAR SPOTS 

HD 189733 is an active K dwarf. Star spots modulate the 
total brightness of the star along the ~12-days rotation cy- 
cle of the star. Coincidentally, the dimming from star spots 
on HD 189733 in the visible, 1-3 %, is comparable to the 
dimming produced by the planetary transit. Accounting for 




Figure 2. Effect of stellar spots on a planetary transit light 
curve. Spots not occulted by the planet produce a deeper transit. 
Spots occulted by the planet produce flux rises on the timescale 
of the transit ingress/egress. 



the presence of star spots when calculating the transmission 
spectrum is therefore essential. Figure [2] illustrates the effect 
of star spots on transit depth measurements. For a plane- 
tary transit across an unspotted star, the depth of the flux 
dimming is proportional to the fraction of the stellar disc 
occulted by the planet, i.e. to the square of the radius ratio: 



d oc 



Api 
A, 



it)' 



If the surface of the star is spotted, then the relation 
becomes: 



d oc 



(1) 



where a is the mean brightness of the part of the star 
occulted by the planet, and a' is the mean brightness of the 
stellar surface compared to a spot-free equivalent. 

Thus there are two distinct ways for star spots to affect 
the recovered radius ratio: 

— first, star spots occulted by the planet during the transit 
reduce the transit depth, leading to an underestimation of 
the transit radius by a factor ^/a. 

— second, star spots not occulted by the planet but sit- 
uated on the side of the star visible during the transit will 
lead to an overestimation of the planetary radius by a factor 



Both effect have a similar amplitude - they are pro- 
portional to the fraction of the stellar flux blocked by star 
spots. Both effects are also wavelength-dependent, because 
star spots have a different temperature and a different spec- 
tral energy distribution than the rest of the stellar disc, so 
that they modify not only the recovered radius ratio, but 
also the inferred transmission spectrum. 

Figure [3] sketches the impact of these two effects on the 
effective transit radius as a function of wavelength (for more 
quantitative details see |Pont et al.|2008| |Sing et al.||2011[ ). 

The present description neglects two further effects of 
variability on the transit light curve: the limb darkening 
of the spots, and the effect of brighter active regions on 
the star, taeniae and plages. The first effect is negligible 
to our level of accuracy. The second could be important, 
but none of the space and ground-based monitoring of the 
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Wavelength 

Figure 3. Sketch of the effect of spots on the inferred transmis- 
sion spectrum, out of a spot occultation (below) and during a 
spot occultation (above). The flat line is the planet transmission 
spectrum inferred in the absence of spots. 



1 US t m En 6 n t 


BJD 


J 




ACS (G800L) 


2453877.20 


-0.0023 


-0.0009 


WFP3 (G141) 


5510.09 


0.0040 


0.0051 


WFP3 (G141) 


5443.52 


0.0024 


0.0074 


NICMOS (F166N) 


4589.36 


-0.0036 


-0.0018 


NICMOS (F187N) 


4571.65 


-0.0050 


-0.0048 


NICMOS (G206) 


4219.98 


-0.0108 


-0.0072 


Spitz (IRACl) 


4429.68 


0.0017 


0.0048 


Spitz (IRACl) 


4039.22 


-0.0159 


-0.0139 


Spitz (IRAC2) 


4427.47 


0.0003 


0.0034 


Spitz (IRAC3) 


4429.68 


0.0017 


0.0048 


Spitz (IRAC4) 


4281.00 


-0.0030 


0.0008 


Spitzer (MIPS) 


4398.60 


0.0047 


0.0063 



Table 3. Values of stellar flux variation used in the spot correc- 
tions, calculated from the Gaussian-process interpolation of the 
APT photometry. 



HD 189733b system available shows a significant signature 
of the crossing of a brighter region. This is coherent with 
our understanding of the variabihty of cool stars and the 
Sun. Dark spots are well-defined, large regions of lower flux, 
while faculae and plages are distributed more evenly across 
the stellar surface and tend to average out in the global flux. 

To recover the planetary radius from the transit data, 
we must estimate a and a' for each transit, including their 
wavelength dependence. 

Fortunately, a large body of data has now been gathered 
on the variability of HD 189733 and the characteristics of 
its spots, and we are in a position to introduce corrections 
for these effects and assess their uncertainty in a much more 
solid way than was possible in earlier studies. 

We describe below how the spot corrections are calcu- 
lated and implemented. It is necessary to delve in some de- 
tails, since they constitute a dominant source of uncertainty 
on the final transmission spectrum. 



3.1 Correction for spots not crossed by the planet 

Spots on the visible surface of the star affect the transmis- 
sion spectrum by causing the occulted and unocculted parts 
of the star to have different spectra. This is the factor a' in 
Equation ^ 

We assume that the mean spectrum of the spots does 
not change with time. Then: 

a' = 1 - f{t)cx 

where f{t) is the flux dimming due to spots at some reference 
wavelength Ao, and cx the ratio of the effect of spot between 
A and the reference wavelength Aq. 



3.1.1 Estimating f{t) 

The factor f{t) is estimated, as in previous studies, by mon- 
itoring the variability of HD 189733. We beneflt from an 
almost continuous photometric monitoring of the star with 
the APT photometer over more than 5 years to measure the 
level of unocculted spots during the HST observations. We 
use Gaussian processes to interpolate the variability data in 
time, as detailed in Appendix A. Gaussian processes allow 



a Bayesian interpolation of the data with a minimum of as- 
sumptions on the functional form and regularity of the light 
curve, and are especially tailored to calculate realistic uncer- 
tainties. The photometric data and interpolation is shown in 
Fig. [i] The resulting factors f{t) at the epochs correspond- 
ing to space observations are given in Table [3] The reference 
wavelength is that used in the APT photometry, the mean 
of the b and y Stromgren filters. 

The ground-based hght curve and its Gaussian-process 
interpolation produce an estimate of f{t) relative to an ar- 
bitrary reference level. To translate that into a value of f{t) 
that can be used to correct the spectroscopic data, it is also 
necessary to know the level of flux corresponding to an en- 
tirely spot-free surface. This reference level cannot be mea- 
sured directly. It does not correspond to the maximum of 
the hght curve, since it is possible - indeed likely - that 
even at its maximum flux the visible surface of the star is 
still affected by spots. 

Fortunately there are several ways to estimate this ref- 
erence level indirectly, and several lines of evidence suggest 
a value of the order of 1-2 % of the flux. 

(i) The statistics of spot-crossing events during the HST 



observations suggests a 1-2 % reference level (see Sing et al 
2011 for details). This is valid if the portion of the star 



crossed by the planet is typical of the whole stellar surface 
in terms of spot coverage. 



(ii) Aigrain et al. (20121 find that in stochastic star spot 



simulations as well as in the Sun, the reference level is above 
the maximum of the hght curve by a value that is compara- 
ble to the variance of the hght curve itself. A much higher 
value requires unnatural spot configurations, and in the case 
of HD 189733 these configurations would have to survive for 
many spot cycles in order to maintain the observed variance 
of the light curve much below the reference level. 

(iii) The transmission spectrum in the visible does not 
show any spectral feature associated with star spots, most 
notably the broad wings of the stellar sodium line and the 
MgH line, which have a steep temperature dependence (see 
Sing et aL]|2011 figure 10). This would be expected if the 
reference level was much higher than 1-2 percent. 
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54000 54500 55000 55500 

Time (MJD) 



Figure 4. Flux measurements for HD 189733 collected with the APT photometer over six years, and Gaussian process interpolation. The 
red line shows the mean value, and the pink area the 1-sigma contour of the possible interpolations. The rectangular shading identifies 
the poorer observing season. The interpolation process is described in Appendix A. 



3.1.2 Estimating c\ 



To estimate C\, we use the same procedure as |Pont et al.| 
( 2008 1 and Sing et al. ( 2011 1. As far as can be gathered from 



the spectrum changes during the spot crossings observed 
with STIS and ACS, a mean temperature of 4250 ± 250 K 
provides a good model of the spectral energy distribution of 
spots. 

We model the spot-free photosphere using a MARCS 
stellar model with solar metallicity, log g = 4.5 and Teff 
= 5000 K. For the spots, we used similar models, but with 
cooler temperatures (4000, 4250 and 4500 K). We refer to 
each spectrum by its temperature, for example F4000(A) is 
the 4000 K spectrum. We compute a model spectrum for a 
star with a fraction S of the projected visible area covered 
in spots at temperature T: 

Fspottedi5,X,T) = (1 - (5)F5ooo(A) + 5FT{\) where we 
have ignored limb-darkening. 

We then compute the resulting flux relative to the spot- 
free case for a given observational setup: 



spotted 



iS,X,T)U{\)dX 



J F5ooo(<5,A,T)L,(A)dA 



where Li (A) is the combined instrument, filter and de- 
tector throughput for observational setup i. Where avail- 
able, we used resolved transmission curves, linearly interpo- 
lated to the resolution of the MARCS spectra. Otherwise 
we used top-hat transmission functions, specifying only the 
minimum and maximum wavelengths. 

To combine the out-of-transit photometry from different 
observatories, we need the ratio between the amplitude of 
the brightness variations in two observational setups. This 
is given by: 



.h{S,T) 



^-.fAS,T) 



We report these ratios in Table [4] where i is always the 
APT. Note that the value of 5 has minimal impact on the 




53500 54000 54500 55000 

Julian Date 



Figure 5. Time series of the difference in flux between the three 
comparison stars in the APT photometry. Note the noisier season 
around JD=55150. 



amplitude ratios, provided 5 « 1. We used 5 — 0.01 to 
obtain the ratios reported in the table. 



3.1.3 Discrepant APT season 

Fig. [5] shows the APT photometry for the three compari- 
son stars relative to each other. This suggests that one of 
the season is of markedly lower quality than the others (be- 
tween JD 2455150 and 2455400) and may have an incorrect 
zero-point for HD 189733. This corresponds to the season 
with the most discrepant behaviour of the Gaussian process 
interpolation for HD 189733 (Fig. [4]| . The Gaussian process 
approach does not include a model of a zero-point change 
(which would be completely degenerate with a change in flux 
of HD 189733). We therefore choose not to use the data for 
this season. 

Only the STIS spectrum was taken during this season, 
and indeed, if we re-scale the STIS and ACS observation 
using the Gaussian process interpolation, we observe a large 
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Instrument 


Wavelength(nm) 


C4000 


C4250 


C4500 


STIS (G430L) 


290 - 


- 370 


1.130 


1.244 


1.403 


RTTS ("043(11^") 


370 - 


- 420 


1.082 


1.160 


1.255 


STIS (G430L) 


420 - 


- 470 


1.049 


1.072 


1.106 


STIS (G430L) 


470 - 


- 520 


1.037 


1.065 


1.087 


STIS (G430L) 


520 ' 


- 570 


0.965 


0.955 


0.946 


STIS (G750M) 


580 - 


- 592 


0.958 


0.930 


0.902 


STIS (G750M) 


592 - 


- 603 


0.936 


0.897 


0.863 


STIS (G750M) 


603 - 


- 615 


0.908 


0.880 


0.856 


STIS (G750M) 


615 - 


- 626 


0.947 


0.908 


0.864 


STIS (G750M) 


626 - 


- 637 


0.917 


0.881 


0.845 


STIS ('G75nM'l 


588 - 


- 590 


1.047 


1.067 


1.049 


STIS fG75nM"l 


589 - 


- 589 


1.086 


1.168 


1.230 


ACS CGsnnLl 


550 - 


- 600 


0.948 


0.922 


0.898 


ACS fGsnnLi 


600 - 


- 650 


0.916 


0.884 


0.850 


ACS fGsnnL'i 


650 - 


- 700 


0.882 


0.849 


0.808 


ACS (G800L) 


700 - 


- 750 


0.811 


0.789 


0.760 


ACS CGSnOLI 


750 - 


- 800 


0.762 


0.746 


0.721 


ACS fGsnnLi 


800 - 


- 850 


0.724 


0.714 


0.691 


ACS fGsnnLi 


850 - 


- 900 


0.702 


0.692 


0.664 


ACS ('GSnnL"! 


900 - 


- 950 


0.675 


0.665 


0.638 


ACS CGSnOLI 


950 - 


- 1000 


0.660 


0.647 


0.619 


ACS fGsnnL) 


1000 - 


- 1050 


0.653 


0.636 


0.603 


WFC3 fG141") 

V V -L V_> y vJX _L rr _L y 


1099 - 


- 1153 


0.601 


0.575 


0.545 


WFC3 I^G141~) 

V V -I- v_/ T J y v_J i rr _L y 


1500 - 


- 1694 


0.484 


0.338 


0.300 


NTCMOS ('C2nfi"l 


1463 - 


- 1520 


0.525 


0.408 


0.365 


NICMOS (G206) 


1520 - 


- 1578 


0.504 


0.365 


0.324 


NICMOS (G206) 


1578 - 


- 1635 


0.470 


0.320 


0.283 


NICMOS (G206) 


1635 - 


- 1693 


0.459 


0.299 


0.263 


NICMOS (G206) 


1693 - 


- 1750 


0.449 


0.290 


0.255 


NICMOS (G206) 


1750 - 


- 1807 


0.447 


0.289 


0.254 


NICMOS (G206) 


1807 - 


- 1865 


0.441 


0.286 


0.250 


NTCMOS (C2n6"l 


1865 - 


- 1922 


0.431 


0.280 


0.244 


NTCMOS ('C2n6'> 


1922 - 


- 1980 


0.447 


0.297 


0.254 


NTCMOS ^0206^ 


1980 - 


- 2037 


0.436 


0.289 


0.250 


NICMOS (G206) 


2037 - 


- 2094 


0.418 


0.276 


0.241 


NTCMOS ^0206^ 


2094 - 


- 2152 


0.411 


0.274 


0.241 


NICMOS (G206) 


2152 - 


- 2209 


0.397 


0.267 


0.235 


NICMOS (G206) 


2209 - 


- 2267 


0.393 


0.270 


0.238 


NICMOS (G206) 


2267 - 


- 2324 


0.388 


0.278 


0.251 


NICMOS (G206) 


2324 - 


- 2381 


0.391 


0.293 


0.267 


NTCMOS (0206^ 


2381 - 


- 2439 


0.392 


0.303 


0.278 


NICMOS (G206) 


2439 - 


- 2496 


0.404 


0.305 


0.272 


NICMOS (F166N) 


1649 - 


- 1667 


0.455 


0.298 


0.263 


NICMOS (F187N) 


1864 - 


- 1884 


0.423 


0.271 


0.236 


Spitzer (3.6) 


3120 - 


- 3190 


0.369 


0.272 


0.231 


Spitzer (4.5) 


4000 - 


- 5000 


0.283 


0.253 


0.230 


Spitzer (5.8) 


4900 - 


- 6300 


0.266 


0.246 


0.224 


Spitzer (8.0) 


6400 - 


- 9300 


0.267 


0.242 


0.210 



Table 4. The scaling factors for the effect of star spots at three 
different temperature, 4000, 4250 and 4500 K. 



discrepancy between the two spectra over the wavelength 
range at which they overlap. Since we do not consider actual 
change in atmospheric properties at this level ("weather") 
to be plausible, we suspect higher uncertainties in the mea- 
surement of the stellar flux during the lower-quality season 
to be responsible for the mismatch. 



Like in Sing et al. (2011 1, we therefore use the assump- 
tion that the spectrum has to be compatible across the over- 
lap in wavelength and choose the spot level that correctly 
connects the STIS and ACS spectra. 



3.2 Corrections for spots crossed by the planet 

Occultations of star spots by the planet are clearly seen in 
the ACS and STIS transit curves. Figures [6] and [7] show the 
residuals around the best-fit transit model (with the same 
system parameters for all data sets), for visible and infrared 
data respectivel^l^ One remarkable feature of the ensemble 
data is that every single visit in the visible shows the sig- 
nature of occulted star spots. No such events are seen in 
the infrared. The signal-to-noise ratio of spot occultations 
is much lower in the infrared data, both because the noise 
on individual data points is higher, and because the spot-to- 
photosphere contrast becomes lower at longer wavelengths. 
This strongly suggests that undetected spot crossings are 
present in the infrared data as well. Figure |8] shows quanti- 
tatively how star spots can be missed in the infrared data. 

At visible wavelengths, we use only the parts of the 
data that appear unaffected by significant star spot cross- 
ings. This corresponds to using a = 1 in Equation [l] 

In the infrared, since individual spot crossings cannot be 
identified, we must rely on an estimate of the average effect 
of occulted spots. If the spots are randomly distributed on 
the surface of the star, the effect of occulted and unocculted 
spots will compensate each other on average (a = a' on 
average). The effect of spot crossing by the planet will be to 
increase the error in this estimate, due to the low-number 
statistics of individual spot crossings. 

We estimate the dispersion introduced by unrecognised 
spot crossings in the infrared in three ways: first using the 
statistics of spot crossings in visible data, second using the 
statistics of repeated depth measurements in the infrared in 
the same passbands, third, using the variability in the light 
curve of the star: 

(i) the standard deviation of the depth measured in the 
ACS and STIS data including the spots (see Figs. |6| is 
314- 10~^. This, while obviously affected by low-number 
statistics, is based on enough transits to be a meaningful 
estimate of the average effect of spot crossings. Scaled with 
the contrasts of Table |4j using T=4250 K for the spots, it 
amounts to a scatter of 94-10"'' on the transit depth near 
1.6 ^m (NICMOS F116N filter) and a scatter of 75-10"^ for 
the Spitzer 8 \im channel. The last number translates into 
0.3 % of the total transit depth. 

(ii) Agol et al. (2008) measured a residual variations of 
0.6 % in the measured transit depth for seven measurements 
at 8 i^m with Spitzer. These authors identify spot crossings 
as a likely dominant source of this scatter, as they note that 
the residuals do not correlate with the stellar flux, indicating 
that the unocculted spots cannot entirely explain this scatter 
(they also argue that random error are smaller than this 
scatter). 

(iii) The total flux of the star in the blue (Stromgren b 
and y filters) varies by 1-3 %, depending on the seasons. 
If the planet does not cross a special latitude respective to 
the position of the spots, then star spots amounting to a 

^ Throughout this paper we designate as "visible" the observa- 
tions shortwards of l|xm, and "infrared" beyond. This does not 
correspond to human eye sensitivity, but rather to the use of 
CCD-detectors for the first category and infrared-array detector 
for the second, which separates the observations into different 
categories for the purpose of reduction and analysis. 
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Figure 6. Residuals compared to the transit model for the STIS and ACS visits, integrated over all wavelengths. The vertical axes are 
scaled so that the same star spot would produce an effect of the same amplitude (using the wavelength scaling in Table |4] with Tspot = 
4250 K). 
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Figure 7. Residuals compared to the transit model for the NICMOS photometric visits, integrated over all wavelengths. The data was 
binned by groups of 7 points to make the temporal coverage comparable to the ACS and STIS data of Fig. [6] The vertical axis are scaled 
so that the same star spot would produce an effect of the same amplitude (using the wavelength scaling in Table|4]with Tspot = 4250 K). 
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Figure 8. Left: residuals around a transit model from the three medium-resolution STIS HST visits, showing the ubiquitous effect of 
occulted spots. Right: Same residuals, with the spot effect scaled to the expected amplitude at 8 i^m, and added errors corresponding 
to those of Spitzer measurements at 8 ^.m (3.7 X 10~* per minute). 



1-3% dimming in the APT b + y passband on the path of 
the planet will scale to a 0.24-0.72 % effect on the depth of 
a transit at 8 [im. 

These three independent estimates are compatible with 
each other, given that the second is the sum of the contribu- 
tion of the occulted stat spots with the effect of unocculted 
star spots (of the same amplitude on average) and the in- 
strumental systematics. A scatter of 0.3 % at 8 \im on a 
due to unrecognised star spot crossings is compatible with 
all three. 

We therefore use 0.3 % at 8 \im as an additional uncer- 
tainty in the depth measurement of individual transits due 
to unidentified spot crossings, and use a =<«'> to calcu- 
late our best estimate of the transit radius in the infrared 
data. This value is scaled with the factors in Table|4]at other 
wavelengths. 



3.3 Limitations of the spot corrections 

Several lines of evidence and cross-checks allow us to build 
some confidence on the spot corrections. We have also stayed 
as conservative as possible in our assumptions about the 
effect of unseen spots in the infrared data. As a results, the 
uncertainties in the infrared measurements are significantly 
larger than in previous studies. 

Nevertheless, some coincidences and compensating ef- 
fects remain possible. We identify some of them here. 

The mean effective temperature of spots is calculated 
from the large spot crossings in the ACS and STIS data. 
It is possible that this temperature is only representative of 
larger spots, and that there is a large population of smaller 
spots with a weaker temperature difference. Such a "leopard 
skin" model for HD 189733 would modify the transmission 
spectrum in ways that would be virtually impossible to cor- 
rect with the available data. Nevertheless, with smaller and 
more numerous spots, the effects of occulted and unocculted 
spots tend to average out over the scale of a full transit. 

The level of unocculted spots, that we have estimated 
to be 1-2% in the previous section, could be much higher. 
This scenario requires a certain number of coincidences to 




Figure 9. Geometry of the HD189733 system. 



remain compatible with the observations. First, the spots 
must be preferentially situated out of the path taken by the 
planet across the star during the transit (the planet crosses 
the 31-55'^ latitude range of the star, see Fig.[9|, so as to re- 
produce the statistics of crossed spots and the variability of 
measured depth at 8 \im. This is possible with a large polar 
spot region, for instance. The spot configuration would also 
have to remain remarkably stable, to reproduce the variabil- 
ity of the APT light curve over more than five years. Also, 
the MgH feature expected with such a high spot dimming 



(see Fig. 10 of Sing et al. 2011 1 is not seen, which requires the 
spot spectrum to be anomalous and different from a cooler 
stellar photosphere (this last point is not as unlikely as it 
sounds, the TrES-1 spot occultation observed with ACS in 



Rabus et al. ( 2009 1 for instance show a spectrum for the spot 



that seems flatter than expected, and in our STIS spectrum, 
the expected MgH spectral feature of the spot near 5000 A 
seems absent). 

The distribution of spots could also be strongly uneven 
in latitude, undermining the a = <a' > relation used for 
the infrared transit. This would bias the connection of the 
infrared and visible parts of the transmission spectrum. 

We consider these situations less likely than our default 
assumptions, and will leave them as possible caveat unless 
they are supported by empirical evidence. 
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3.4 "Weather" variability 

Reconstructing a unique transmission spectrum from data 
taken at different times rests on the assumption that in- 
trinsic variations in the planetary spectrum ( "weather" ) are 
smaller than the observational uncertainties. 

Any variation between different runs at the same wave- 
length could be explained by actual variations in the trans- 
mission spectrum of the atmosphere of the planet. Neverthe- 
less, in our most accurate datasets (ACS, STIS, NICMOS 
filters), radius measurements repeated with the same instru- 
ment agree with each other within the uncertainties without 
intrinsic variations in the planet (in the infrared, the uncer- 
tainties related to spot and instrument effect corrections are 
larger than the variations between different measurements). 

Moreover, the radius variations that we measure span 
several atmospheric scale heights, which is much larger than 
the expected variations due to planetary weather. 

We therefore assume that the observed variations are 
due to the mean atmospheric transmission spectrum around 
the planet combined with observational uncertainties, rather 
than the actual variations of the atmospheric transmission 
spectrum. 

Building a single transmission spectrum from our data 
sets rests on this assumption. 



4 RESULTS 

4.1 Combined UV-to-IR transmission spectrum 

Table[5]and Figure [To] give the resulting values for the planet 
size as a function of wavelength, after applying the spot 
correction and uncertainty estimates detailed in the previous 
Section, using Tspots = 4250 ± 250 K. 

In Table[6]and Figure [TTj the spectrum data is averaged 
over a limited number of passbands, arbitrarily chosen when 
combining data from different instruments. A condensate- 
free model of the transmission spectrum of HD 189733b from 



Passband 
(nm) 



Radius ratio 



Fortney et al. ( 2010| is plotted for comparison. 

Shortward of 1 /.im, the broad transmission spectrum 
is very well measured, and defined by a single steep blue- 
ward slope. The core of the sodium doublet is resolved by 
the higher-resolution STIS data, that also shows an absence 
of broad wings. The implications were already presented in 



Pont et al.N2008|),|Lecavelier Pes Etangs et al.| ( |2008| ), |Sing 
et al.|(|201l| andlHuitson et al.l (|2012l). The overall slope is 



compatible with Rayleigh scattering by solid or liquid par- 
ticles (with sizes below 0.1 (im), with a A"** dependence 
of the cross section, evenly distributed in an atmosphere 
at T ~ 1300K. The increasing slope towards the UV and 
the height of the sodium line suggest a temperature rise 
in the upper atmosphere. The core of the sodium doublet, 
and possibly the sodium doublet, are visible above the haze 
Rayleigh signature and constrain the altitude at which the 
haze is seen: low enough to leave the core of the lines seen, 
but high enough to cover the pressure-broadened wings of 
the sodium and potassium doublets. 

The infrared data is compatible with a featureless spec- 
trum, as well as with the presence of muted molecular fea- 
tures. With the addition of the uncertainties for undetected 
spot crossings, none of the measured variations in radius 
in the infrared exceeds the uncertainties. We note that our 



320-370 





15811 







00032 


370-420 





15751 


± 





00016 


420-470 





15718 


± 





00014 


470-520 





15695 


± 





00013 


520-560 





15666 


± 





00013 


560-580 





15644 


± 





00014 


580-595 





15679 


± 





00007 


595-600 





15631 


± 





00022 


600-650 





15609 







00009 


650-700 





15585 


± 





00011 


700-750 





15572 


± 





00011 


750-800 





15586 







00012 


800-850 





15552 







00012 


850-900 





15553 


± 





00013 


900-950 





15546 


± 





00013 


950-1000 





15552 


± 





00016 


1000-1170 





15512 


± 





00022 


1450-1750 





15480 


± 





00026 


1750-2100 





15491 


± 





00025 


2100-2500 





15461 


± 





00030 


3200-3900 





15507 


± 





00027 


4000-5000 





15542 


± 





00035 


5000-6400 





15476 


± 





00067 


6400-9400 





15510 


± 





00034 


23500-24500 





15459 


± 





00094 



Table 6. Radius ratio as a function of wavelength, with spot 
effects accounted for, and binned in 25 wavelength intervals. 



multi-instrument, multi-transits approach is not best suited 
to the detection of features in specific spectral ranges. This 
is better done by analysing data from a single instrument, 
preferably acquired during a single transit. The issue of in- 



frared features is discussed in detail in Gibson et al. (2011 1 
for the 1-2 iJ.m range and |Desert et al.| ( |2011[ ) at longer wave- 
lengths. 

It is remarkable, though, that the more data sets ac- 
cumulate, the nearer they evolve towards a continuation of 
the featureless slope indicated in the visible. Indeed, the 
most economical inference from the observations would be 
a monotonic decrease of the transit radius towards longer 
wavelengths, with no spectral feature rising above the noise 
level. 

By contrast, clear-atmosphere models struggle to repro- 
duce the observed transmission spectrum of HD 189733b. 
This is most obvious in the visible, but is also the case in 
the infrared. Clear-atmosphere models with molecules pre- 
dict water features between 1 and 2 [xm that are incompat- 



ible with the NICMOS filter data ( |Sing et al.||2009) , and 
a rising opacity and larger transit radius at 8 [xm than at 
4.5 \i.ra, that the data do not suggest. An entirely flat spec- 
trum in the infrared produces a lower to the data than 
clear-atmosphere models. 



4.2 Simple extended-haze model 

The 0.6-1 \im transmission spectrum from ACS was inter- 



preted by Lecavelier Des Etangs et al. ( 2008 1 and subse- 



quently as the signature of a single high-altitude layer of 
scattering haze. The grains must be abundant enough to 
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Instrument 


BJD 


Central 


Band 


Rp/ Rs 




(setting)visit 


[-2450000] 


Wavelength [nm] 


half-width [nm] 










STIS(G430L)1 


5155.12490 


3300 


400 


0. 


15866 





.00043 


STIS{G430L)2 


5334.80977 


3300 


400 


0. 


15734 





.00051 


STIS(G430L)1 


5155.12490 


3950 


250 


0. 


15762 





00020 


STIS(G430L)2 


5334.80977 


3950 


250 


0. 


15732 





00026 


STIS(G430L)1 


5155.12490 


4450 


250 


0. 


15701 





.00024 


STIS{G430L)2 


5334.80977 


4450 


250 


0. 


15728 





.00018 


STIS{G430L)1 


5155.12490 


4950 


250 


0. 


15669 





00023 


STIS(G430L)2 


5334.80977 


4950 


250 


0. 


15706 





.00015 


STIS(G430L)1 


5155.12490 


5450 


250 


0. 


15654 





00023 


STIS{G430L)2 


5334.80977 


5450 


250 


0. 


15672 





.00016 


ACS(G800L)1 


3877.20896 


5750 


250 


0. 


15644 





.00014 


STIS{G750M) 


5148.46821 


5865 


56 


0. 


15638 





00027 


STIS(G750M) 


5148.46821 


5895 


11 


0. 


15703 





.00011 


STIS(G750M) 


5148.46821 


5980 


57 


0. 


15631 





00022 

,\J\J\J ^ ^ 


STIS{G750M) 


5148.46821 


6095 


55 


0. 


15617 





00036 


STIS{G750M) 


5148.46821 


6207 


57 


0. 


15600 





00027 


STIS(G750M) 


5148.46821 


6321 


57 


0. 


15611 





.00019 


ACS(G800L)1 


3877.20896 


6250 


250 


0. 


15610 





00012 


ACS(G800L)1 


3877.20896 


6750 


250 


0. 


15585 





.00011 


ACS(G800L)1 


3877.20896 


7250 


250 


0. 


15572 





00011 


ACS(G800L)1 


3877.20896 


7750 


250 


0. 


15586 





.00012 


ACS(G800L)1 


3877.20896 


8250 


250 


0. 


15552 





00012 


ACRCG8nnT.'n 


3877 20896 


8750 


250 


0. 


15553 





00013 


ACRCG8nnT.'n 


3877 20896 


9250 


250 


0. 


15546 





00013 


ACS(G800L)1 


3877.20896 


9750 


250 


0. 


15552 





.00016 


ACS(G800L)1 


3877.20896 


10250 


250 


0. 


15496 





00024 


HSTfWFC'il2 


5510 09830 


11050 


230 


0. 


15671 





00084 


HSTfWFC^'ll 


544S 52370 


11335 


345 


0. 


15549 







NICMOS(G206) 


4219 98178 


14920 


285 


0. 


15335 





00129 


NICMOS(G206) 


4219 98178 


15500 


285 


0. 


15330 





. uuuoo 


HSTfWFC'^'12 


5510 09830 


15925 


785 


0. 


15608 





00055 


HSTfWFC^'ll 


544S 52370 


16070 


860 


0. 


15543 





. uuuo^ 


NICMOS(G206) 


4219 98178 


16070 


285 


0. 


15241 





.00078 


NICMOS(F166N)l 


4589 36894 


16580 


90 


0. 


15516 







NICMOS(F166N)2 


461 1 60656 


16580 


90 


0. 


15582 







NICMOS(G206) 


4219 98178 


16650 


285 


0. 


15405 





.00079 


NICMOS(G206) 


4219 98178 

rr _L C • 1.^ O J- 1 


17220 


285 


0. 


15424 







NICMOS(G206) 


4219 98178 


17790 


285 


0. 


15511 





. uuuuo 


NICMOS(G206) 


4219 98178 


18370 


285 


0. 


15490 





.00070 


NICMOS(F187N)l 


4571 65062 


18740 


100 


0. 


15459 





. uuuou 


NICMOS(F187N)2 


4689 26721 


18740 


100 


0. 


15456 





.00048 


NICMOS(G206) 


4219 98178 


18940 


285 


0. 


15572 





.00075 


NICMOS(G206) 


4219 98178 


19510 


285 


0. 


15525 





. \j\J\J\JZj 


NICMOS(G206) 


4219 98178 

^ _L C • 1^ O J- 1 O 


20090 


285 


0. 


15370 





. uuuoo 


NICMOS(G206) 


4219 98178 


20660 


285 


0. 


15486 





00080 


NICMOS(G206) 


4219 98178 


21240 


285 


0. 


15395 





nnnsn 


NICMOS(G206) 


4219 98178 


21810 


285 


0. 


15483 







NICMOS(G206) 


4219 98178 


22380 


285 


0. 


15490 





.00061 


NICMOS(G206) 


4219 98178 


22960 


285 


0. 


15491 





.00073 


NICMOS(G206) 


4219.98178 


23530 


285 


0. 


15432 







NICMOS(G206) 


4219.98178 


24110 


285 


0. 


15496 





.00070 


NICMOS(G206) 


4219.98178 


24680 


285 


0. 


15520 





00090 


Spitz{IRACl) 


4429.68978 


36000 


4000 


0. 


15471 


0, 


.00051 


Spitz(IRACl) 


4039.22278 


36000 


4000 


0. 


15547 


0, 


.00037 


Spitz{IRACl) 


5559.55455 


36000 


4000 


0. 


,15452 


0, 


.00059 


Spitz{IRAC2) 


4427.47301 


45000 


5000 


0. 


,15538 


0, 


.00051 


Spitz(IRAC2) 


5189.05249 


45000 


5000 


0. 


15543 


0, 


.00049 


Spitz(IRAC3) 


4429.68978 


58000 


7000 


0. 


15476 


0, 


.00067 


Spitz{IRAC4) 


4281.00701 


78500 


14500 


0. 


15510 


0, 


.00034 


Spitzer(MIPS) 


4398.60560 


240000 


45500 


0. 


,15459 


0, 


.00094 



Table 5. Transmission spectroscopy results for all data sets. Note that the uncertainties on the transit radius arc not independent. 
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Figure 10. Transmission spectrum data, with data sets and visits indicated separately. Lines connect values obtained at the same time 
with the same instrument. 
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Figure 11. Our combination of the available data into a single set of constraints, assuming that the intrinsic transmission spectrum of 
the planet is constant with time without our uncertainties. The grey line shows a synthetic spectrum with a dust-free model. 



provide a higher opacity than the wings of the sodium and 
potassium lines, and transparent enough in the visible so 
that the Rayleigh slope dominates over absorption. Dust 
and clouds have been considered as a possible important 
components in the atmosphere of hot Jupiters since the first 
atmosphere models and observations. The atmospheric tem- 
perature of hot Jupiters correspond to spectral type L for 
brown dwarfs, and observations have shown that these ob- 
jects have very red colours, interpreted as due to the effect 



of dust in their atmosphere. Several common grain- forming 
elements have condensation temperatures above the tem- 
perature of late L-type objects and most hot Jupiters, in- 
cluding enstatite (MgSiOs), forsterite (Mg2Si04), corundum 
(AI2O3) and elemental iron (Fe). Of these, only enstatite is 
transparent over the visible wavelength range. 

The explanation of the ACS result in terms of a thin, 
high-altitude layer of enstatite haze made specific predic- 
tions about the rest of the spectrum that were not borne 
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out by the observations. The spectrum was expected to flat- 
ten towards the UV, as the wavelength becomes comparable 
to the size of the grains and therefore moves to the flat part 
of the Mie scattering curve. It was also expected to drop 
sharply towards the infrared, where the line of sight reaches 
regions below the haze layer. Measuring these two features 
was the main motivation of the HST observing programme 
GO-11740. Neither was corroborated by the new observa- 
tions with STIS and WFC3. 

The ACS and STIS observations are compatible with 
the dependence of Rayleigh scattering over more than 
5 atmospheric scale heights, i.e. more than two orders of 
magnitude in pressure (see Figure . In the infrared the 
spectrum becomes flatter, and remains above the extension 
of the Rayleigh slope. 

An extension of the haze hypothesis can account for 
these observations. Condensates not conflned to a high- 
altitude layer, but extending over most of the atmosphere, 
can be expected to produce this type of spectrum. At higher 
altitudes, small grains produce the Rayleigh slope. At lower 
altitudes, grains become larger, and their size gets closer to 
the wavelength, producing a flatter spectrum, until a cloud 
deck layer is reached that entirely scatters or absorbs the 
incoming star light. 

A full model of the transmission spectrum of HD 
189733b in terms of haze and clouds is beyond the scope 
of this paper, and probably unwarranted given the present 
amount of observational constraints. We develop here a flrst- 
order, qualitative model to account for the observed spec- 
trum with few assumptions and free parameters. 



Helling (20031, we estimate the settling timescale by calcu- 



lating the time taken for grains to cross one atmospheric 
scale height at their terminal fall velocity.The terminal ve- 
locity of the particles is obtained by balancing the pull of 
gravity against the drag of the gas flow: 



is the atmospheric scale height (fc is the 



where H = - 
/J 

Bolzmann constant, T the temperature, ^ the mean molec- 
ular weight, rriu the atomic mass unit, g the gravity in the 
atmosphere), iifaii is the terminal velocity for particles of a 
given size due to atmospheric drag. 

Different expressions for the atmospheric drag can be 
used depending of the flow regime affecting the grains. For 
low gas densities, the flow is molecular (Brownian motion). 
The transition from molecular to viscous flow is described 
by the Knudsen number, deflned as the ratio of the mean 
free path to the particle size: 



I 

2a 



Values of Kn near 1 separate molecular flow, where indi- 
vidual impacts dominate the dynamics of the particle, from 
laminar flow, where the gas is dense enough to be treated 
as a viscous fluid. 

For the grain sizes, gas pressures and temperatures rel- 
evant here (a < 10 (xm, p < 1 bar,r 1000 K), the Knudsen 
number is much larger than unity, so that the flow is closer 



to molecular than to viscous. In that case ( Woitke & Helling 



2003h: 



4-2.1 Rayleigh slope 

The Rayleigh scattering cross section evolves as: 



a 
A* 



where a is the grain size and n the refraction index. Since 
the dependence on grain size is very steep , a oc a®, for most 
grain size distributions, the largest grains will dominate the 
scattering cross section. In the case of HD 1897833b, the 
transmission spectrum implies that the maximum grain size 
remains similar over five scale heights or more in the upper 
atmosphere. It also implies that their abundance remains ap- 
proximately constant. Therefore, keeping the Rayleigh slope 
in the transmission spectrum over such a large wavelength 
range requires well-mixed grains with a constant size cutoff. 



4-2.2 Settling regime 

The size distribution of grains with altitude for condensates 
in a planetary atmosphere is controlled by a balance between 
the timescale for the gravitational settling of the grains, and 
the timescale for replenishment of the grain population. 

At a given altitude, the maximum grain size a^ax sus- 
tained will be the one for which the two timescales coincide: 

'^sct (f^max) — '^rcp(ctmax) 

where Tsct and Trcp are the settling and replenishment 
timescales. 



Following Ackerman & Marley (20011 and Woitke & 



I TT ^ pcond g 
4 Pp-as Ct 



where pcond and pgas are the density of the grains and gas, 
and Ct is the sound speed in the gas: 



Cf = 



-ykT 
jj.mu 



Usi ng the ideal gas equation of state, Pcond = 3.2 x 10' 



kgm (Ackerman & Marley 



2001 



for enstatite) and the 
following parameters for HD 189733b, T=1200 K, g = 21 
ms-^ ^1 = 2.35, 7 = 7/3, gives for molecular flow: 



Tsct 



P 



[Ibar] \ [l^m] 



10*s 



Thus, for instance, the settling timescale for 0.1 [im at 
10 mbar is about one month. 

For the purpose of a flrst-order understanding of the 
transmission spectrum, the important parts of the expres- 
sions above for wfaii is its dependence on particle size and 
gas pressure: 



"fall 



1 -1 

^ Pgas 



The replenishment timescale Trop cannot be estimated 
robustly from present models and observations. Equilibrium 
1-D atmosphere models do not predict the formation of 
any condensates, because the atmosphere is vertically stable 
against convection. 3-D models tracking the motion of grains 
have not been constructed yet for hot Jupiters, although this 
is expected to change soon. 
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Two processes could account for grain formation: a ver- 
tical exchange of mass in the atmosphere, due to the vigor- 
ous day-night atmospheric circulation, for instance via eddy 
diffusion, or the formation of photochemical haze through 
the direct action of the stellar irradiation. In both cases, 
the grain formation timescale is not known at present. We 
therefore leave Trop as a free parameter. 

In the case of vertical mixing, note that this parameter 
can be directly linked to the vertical momentum diffusion 
parameter Tfzz, used for instance in |Spiegel et al.| ( |2009[ ) and 



Youdin & Mitchell ( 2010 1 to parameterise vertical mixing in 



hot Jupiters: 

- I!l 



dwarf models ( [Freytag et al.|2010[ ). [Spiegel et al.| ( |2009[ ) find 
that values up to 10 cm s~ are required to keep TiO aloft 
in hot Jupiters and produce a stratospheric inversion. 

A value of K^z ~ 10^^ corresponds to a replenishment 
timescale of 1 hour for HD 189733b, a value of K^z = 10^ to 
10* years. Note that there is no physical reason for K^^ to 
be constant as a function of height in the atmosphere. On 
the contrary, large changes are to be expected. Convective 
mixing for instance tends to produce Kzz oc p~^^^. 

If the formation process produces grains of all sizes, and 
the size distribution at any given height is dominated by 
the balance between settling and replenishment, then in the 
molecular regime, Tset oc ap~^ implies, assuming a constant 

amax(p) OC p 

This relation between particle size and pressure corre- 
sponds to: 

d7? = -4/7dln(A) 

if the initial grain size distribution is assumed to be flat in 
terms of number of grains (n(a) oc 1), and 

d7? = -dln(A) 

if the initial grain size distribution is assumed to be flat in 
terms of mass fraction of grains of difference sizes (n(a) oc 
The derivation of these relations is given in Ap- 



-1/3N 



pendix C. 



4.2.3 Cloud deck 



A third possible regime is a "cloud deck" dominated by large 
grains, which produces a linear cut at a given height in the 
transmission spectrum. 



4.2.4 Haze+cloud scenario 

Figure [12] illustrates these three regimes in the transmission 
spectrum, (i) Rayleigh scattering by well-mixed small grains, 
(ii) settling grains, and (iii) cloud deck. It illustrates how a 
haze-|-cloud scenario could account for the broad features of 
the observed spectrum. 

This haze-|-cloud scenario is partly inspired from exam- 
ples in the Solar System. Venus and Saturn for instance have 
an atmosphere dominated by clouds, with a layer of haze 



above the clouds made of smaller particles, that dominates 
the opacity to incoming sunlight in the visible. 

We note that the observed transmission spectrum sug- 
gests a "bottleneck" effect on K^^. From a certain point 
upwards, the vertical mixing must grow faster than the set- 
tling timescale (thus faster than p^^) to keep the maximum 
grain size constant. In that case, the largest grains that can 
stay aloft above the cloud deck (~ 0.1 \x.m in size) are well 
mixed all the way to the highest layers. It seems to be the 
only way to produce a constant Rayleigh slope over several 
scale heights, as observed. 



5 HAZE AND CLOUDS IN THE PHASE 
CURVE AND ECLIPSE DATA 

5.1 New/revised eclipse and phase curve data 

Information on the atmosphere of HD 189733b is also avail- 
able from measurements of the secondary eclipse and phase 
curve in several infrared channels. The secondary eclipse in 
a given passband measures the brightness temperature in- 
tegrated over the whole day hemisphere of the planet. The 
phase curve measures the day-night contrast in brightness 
temperature. The phase curve also measure the shift of the 
hottest spot on the planet, due to the jets redistributing 
heat from the day side to the night side. 

The realisation that the transmission spectrum from 
UV to mid-infrared is very different from clear-atmosphere 
predictions leads us to revisit these eclipse and phase curve 
data as well. It would be surprising that the haze/cloud/dust 
component were to dominate the transmission spectrum en- 
tirely around the whole planetary limb, and have no effect 
on the emission spectra measured by the eclipse and phase 
curve. Condensates will not only affect the shape of the 
emergent spectrum, but by modifying the radiative trans- 
fer, affect the redistribution of heat around the planet as 
well. 

hereafter K12) presents new 



Knutson et al. (2012 



eclipse and phase curve data, and a re-discussion of previ- 
ous data sets, in the Spitzer passbands at 3.6, 4.5, 8 and 24 
(im. Their results are summarised in Table [7| The seminal 
results of Knutson et al. ( 2007b I showed that at 8 \xm, the 



Spitzer light curve indicates an eastward shift of the hot 
point of 16 ± 6°, and a day-night temperature contrast of 
~350 K (K12 show that instrumental systematics make this 
last number unreliable). The phase curves at 3.6, 4.5 and 
24 (j.m present eastward shifts (listed in Table [7| that at 
compatible with the finding at 8 \i.m. 



5.2 Puzzling features explained by the 
haze/clouds scenario 

Intriguingly, there are several features of the secondary 
eclipse and phase curve information that the clear- 
atmosphere models in K12 struggle to explain, and that 
would be natural consequences of an atmosphere dominated 
by haze opacity. 

K12 do not include the condensates in their modelling 
and discussion of the eclipse and phase curve data. The 
reason, presumably, is that the emission spectrum probes 
deeper regions of the atmosphere than the transmission 
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Figure 12. Spectrum and model as in previous figure. The dotted lines, from left to right, indicate the effect of Rayleigh scattering at 
2000K, 1300 K, a cloud with grain sizes increasing linearly with pressure, and an opaque cloud deck. 



Wavelength 


AT 


^day 


^hotspot 


'/'hotspot 


[im 


[K] 


[K] 


observed 


model 


3.6 


503 (21*) 


1328 (11) 


35° (6) 


53° 


4.5 


264 (24*) 


1192 (9) 


20° (5) 


52° 


8 


350 (:) 


1259 (7) 


23° (3) 


49° 


24 


188 (48) 


1202 (46) 


37° (8) 


49° 



Table 7. Day-night temperature contrast, day-side temperature, 
and hot-spot eastward longitude shift, observed and expected 
from a circulation model, for HD 189733b. The temperatures 
are geometric means over the planetary disc. Data from |Knut-| 
[son et ah] (|2012|l except th e 8 \im and 24 |j.m day-side tempera- 
tures from |Cowan &: Agol| l [2011[ |. (*) according to our re-analysis 
in Appendix B, the uncertainties of the day-night temperature 
differences at 3.6 and 4.5 |im should be it 150 K and ±78 K. 



spectrum (by a factor ~ 50, because of the grazing geome- 
try), and that they assume that the haze layer forming the 
transmission spectrum in the visible is confined to a high 
altitude and transparent in the infrared. However, the trans- 
mission spectrum data presented in this paper suggests that 
the haze extends over more than a factor 100 in pressure, 
and that it may be a source of opacity in the infrared as 
well. 

The first-order effect of dust/haze/clouds is to move 
the contribution functions and effective photosphere at all 
wavelengths to higher altitudes (lower pressure). Based on 
the transmission spectrum, we can expect this to be stronger 
at shorter wavelengths. Haze does not suppress atmospheric 
features in the emission spectrum entirely, but reduces their 
amplitude, because scattering increases the optical path for 
a given pressure difference. This is well studied in the spec- 
trum of L-type brown dwarfs. By moving the photosphere 



to lower pressures, haze and clouds will also tend to re- 
duce the timescale of heat loss from the atmosphere. In a 
hot Jupiter, this will tend to decrease the efficiency of the 
dayside-to-nightside heat redistribution, which can be ob- 
served as an increased day-night temperature contrast and 
decreased eastward drift of the hottest point in the atmo- 
sphere. 

Four possible signatures of the presence of dust in the 
eclipse and phase curve information for HD 189733b are 
listed below. A comparison of the capacity of the "clear" 
and "dusty" scenarios to account for the observations is sum- 
marised in Table IHl 



5.2.1 Larger phase curve amplitude at 3.6 microns than 
in other passbands 

Atmosphere models predict that the amplitude of the day- 
night temperature contrast decreases with altitude, because 
the radiative timescale increases with higher densities much 



faster than the speed of the wind decrease (see e.g. Showman 
fc Guillot||2002| [Cooper fc Showman||2005[ ). In such case, a 
larger day-night temperature contrast indicates that we see 
the flux coming from a higher layer in the atmosphere. 

The day-night temperature contrast at 3.6 (im is al- 
most double that measured in other Spitzer passbands (see 
Table [7|. K12 note that this is puzzling for a clear atmo- 
sphere dominated by absorption from molecular bands. In 
any non-pathological configuration, a higher temperature 
contrast indicates a higher-altitude contribution function. 
But all models based on molecules abundant in a hot Jupiter 
atmosphere suggest that the opacity at 4.5 |i.m should be 
higher than at 3.6 \im. The 4.5 [im band should therefore 
sample a higher layer of the atmosphere, and show a larger 
temperature contrast. The opposite is observed. 
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This is arguably the observation in the Spitzer data that 
is the most robust and tightly connected to the physics of 
the atmosphere. It is unexplained in K12. It also affects the 
interpretation of the day-side emission spectrum, because 
the assumption of a higher opacity at 4.5 ^.m is the founda- 
tion of the measurement of the temperature inversion in the 
atmosphere from the Spitzer passbands. 

However, the presence of condensates can make the at- 
mosphere more opaque at 3.6 \im than at 4.5 [im in spite 
of molecular absorption. A decrease of opacity with wave- 
length extending to the near infrared is suggested by the 
transmission spectrum. 



5.2.2 Reduced hot spot shift 

The longitude of the eastward shift of the hot spot is ex- 
pected to result from a competition between the heat trans- 
fer timescale and the heat loss timescale (identified in hot 
Jupiters to the advective timescale Tadv and the radiative 
timescale Trad, e.g. Showman & Guillot 2002). 

The observed eastward shift of the hot spot in the at- 
mosphere is about half the expected values in the Spitzer 
channel. K12 discuss several possible explanations: a planet 
rotating slower than tidal circularisation, enhanced heavy- 
element abundances in the atmosphere, and an increased 
atmospheric drag. A slow-rotating planet or increased drag 
would increase Tadv by slowing the atmospheric heat redistri- 
bution, whereas more heavy elements would reduce Trad by 
allowing more rapid cooling through atomic line emission. 

An extended haze cover, however, would naturally bring 
the hot spot nearer to the sub-stellar point. Increased opac- 
ities move the photosphere to lower pressures, where the ra- 
diative timescale is shorter. To first order the eastward drift 
of the hot spot depends on the ratio Trad/Tadv Trad changes 
much more rapidly than Tadv with pressure in hot Jupiters, 
to first order Trad ~ in the upper atmosphere, whereas 

Tadv ^ CSte. 

Therefore, raising the photosphere by about one scale 
height will result in decreasing the hot spot shift by 1/e. 
This is compatible with the indications from the transmis- 
sion spectrum. The infrared opacities are roughly one scale 
height higher in the transmission spectrum than the clear- 
atmosphere expectations, which would predict a latitude for 
the hot spot reduced by the same factor compared to the 
clear-atmosphere models, as suggested by the observations. 



5.2.3 Large brightness temperatures at 8 \im and 24 \im 

K12 conclude that the observed 8 iJ.m and 24 \im brightness 
temperatures are larger than the values predicted by clear- 
atmosphere models, both on the day side and on the night 
side. They are nearer to a black-body curve than to the 
expected spectrum sculpted by molecular absorption lines 
(see Fig. 



13 1. K12 find the observations impossible to fit 



even with the added free parameters of a high-altitude grey 
absorber and non-equilibrium chemistry. Only modifying the 
planet rotation rate and the elemental abundances produces 
a satisfactory fit. 

However, a brightness temperature at 8 and 24 \im 
closer to the shorter wavelengths than for a clear-atmosphere 
models is a natural prediction of a haze and cloud-deck 



configuration. Scattering and absorption by condensates 
shorten the altitude differences sampled by the different 
wavelengths. For a given temperature gradient with altitude, 
the brightness temperatures will be correspondingly closer, 
bringing the emergent spectrum nearer to a blackbody. 



5.2.4 Flat spectrum at 5-14 [i-m 

The observations at 3.6, 4.5, 5.8, 8, 16 and 24 \im with 
Spitzer are passband-integrated time series, and therefore 
do not provide information on the local shape of the emis- 
sions spectrum. Grillmair et al.| ( |2007| ) observed the eclipse 
spectrum between 5 and 14 \im with the IRS spectrograph 
on Spitzer. These results show a less marked feature than in 
clear-atmosphere models (see figure 10 of K12, repeated in 



Fig. 13 1. Molecular features are predicted over this interval 



by clear-atmosphere models, but are found to be muted in 
the data. In best-fit model from K12, for instance, the uncer- 
tainties must be stretched to make the observed spectrum 
marginally compatible with the model. While this configu- 
ration produces a satisfactory statistics, it is improbable: 
random errors do not usually conspire to erase spectral fea- 
tures. 

Again, lower-amplitude features are a natural predic- 
tion of cloudy/haze models. Reducing the amplitude of 
molecular feature is an observed effect of condensates in L- 
type brown dwarf spectra. 



5.3 "Dusty" scenario for HD 189733b 

Altogether, the simplest scenario compatible with the trans- 
mission spectrum - i.e. a gentle overall decrease of the opac- 
ity with wavelength - also explains observed anomalies in the 
emission and phase curve data. A global haze/cloud, with 
grains large enough near the infrared photosphere to affect 
radiative transfer, moves the photosphere upwards, reduces 
the hot spot shift, and damps the differences between the 
depth of the contribution functions for the different Spitzer 
channels. It may also make the atmosphere more opaque at 
3.6 (im than at longer wavelengths. 

K12 explain the above features in a piecemeal fashion, 
invoking several additional ingredients, such as enhanced 
metallicity, non-equilibrium chemistry, absence of tidal syn- 
chronisation. All these explanations are plausible. Still, none 
of the clear-atmosphere models can explain the high day- 
night contrast at 3.6 ^m, which appears to be one of the 
most robust Spitzer observations. 

The introduction of "dusty" models was a leap forward 
in the study of brown dwarfs, and our results suggests that a 
similar effort may be useful for hot planets like HD 189733b. 



5.4 Remaining discrepancies 

Some features of the observations do not fit neatly in the 
simplest models, even when including the possible effects of 
condensates. 

One is the marginal indication in the transmission spec- 
trum that the radius at 4.5 ^.m is larger than at 3.6 \i.m. 
Uncertainties due to the effect of star spots on the transit 
depth (discussed in Section |3| imply that this is only at the 
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Observation 


Dusty atmosphere 


Clear atmosphere 


Shortwave transmission spectrum 


ok 


no 


Higher day-night contrast at 3.6 |J^m 


ok? 


no 


Longwave day-side brightness 


ok 


no 


Featureless IRS spectrum 


ok 


difficult 


Low hot-spot shift 


ok 


difficult 


High transmission opacity at 4.5 |J.m 


difficult 


ok?* 



Table 8. Compared capacity of "dusty" and "clear" atmosphere scenarios to account for observations. *: requires non-equilibrium CO 
absorption according to K12. 



1 — a significance level, although both Desert et al. ( 2011 1 



and K12 observe a larger radius at 4.5 |xm. 

This would be difficult to reconcile with the phase curve 
amplitude being larger at 3.6 \im. Only large changes in 
the relative opacities between the two channels around the 
planet could account for both features. 

There is also a marginal indication of a lower opacity 
at 3.6 \im from the higher hot spot shift. The simplest haze 
model would predict the phase shift to decrease with in- 
creased opacity. This is a less stringent constraint than the 
day-night temperature contrast though, because the lon- 
gitude of the hot spot depends on the details of the at- 
mospheric circulation, and the Tadv/Tiad approximation is 
rather crude. 

We note that there is a tension between the phase curve, 
hot point shift and transit depth results at 3.6 and 4.5 \unm 
K12, in both the clear and dusty scenario. Given the impor- 
tance of star spots and instrumental effects, it is possible 
that one of these results will turn out to be modified in light 
of further observations. 

K12 give the highest significance to the phase curve 
amplitude, and consequently this observation has a strong 
weight in our interpretation. However, the phase curve of 
the planet is entangled with the variability of the hot star 
in the infrared light curves. K12 use the Gaussian-process 
interpolation of the APT light curve given in Appendix A to 
constrain the variation in the stellar flux during the plane- 
tary orbit to recover the phase curve. There are two pitfalls 
with this procedure though. The first is that there is no APT 
measurements close in time to the Spitzer measurements, so 
that a large extrapolation over time is required. The second 
is that the time sampling of the APT light curve (one point 
per day typically) is not adapted to resolving the curvature 
of the light curve on the time scale of the planetary orbit 
(2.2 days). This implies that the constraint on the shape of 
the stellar contamination to the measured phase curve will 
depend mainly on the assumptions made on the covariance 
kernel of the Gaussian process, rather than on the data itself. 

In Appendix B, we study how these two factors may 
have led to an under-estimation of the uncertainties on the 
phase curve amplitudes at 3.6 and 4.5 [im in K12. We are 
helped by the availability of one month of precise continuous 
monitoring of HD 189733b with the MOST sateUite. This 
provides an empirical check on the capacity of the Gaussian- 
process interpolation to constrain the curvature of the stellar 
light curve over the time scale of the planetary orbit. 

We find that indeed, the capacity of the Gaussian- 
process interpolation to constrain the stellar contribution is 
severely limited. Consequently, we recalculate the uncertain- 



ties on the phase curve amplitude as 27-10"'' and 25-10"^ for 
3.6 and 4.5 [im respectively (instead of 6- 10~^ and 9-10"^ in 
K12). This translates in uncertainties of ±150 K and ±78 K 
on the temperature contrasts (instead of ±21 K and ±24 K). 
The details are given in Appendix B. With these uncertain- 
ties, the phase curve amplitude at 3.6 \im is still larger than 
at 4.5 [i.m, but with a significance reduced to ~ 2a. 



6 A NEW PICTURE OF THE ATMOSPHERE 
OF HD189733B 

In this Section we explore some consequences of the 
haze±cloud scenario for the planet HD189733b. 



6.1 Is a hot stratosphere ruled out? 

The presence of an inverted temperature profile near the 
photosphere has been inferred for several hot Jupiters from 
the relative fluxes in the 3.6 and 4.5 \im Spitzer bands ( |Fort- 
ney et al. 20101. Since the opacity at 4.5 [im is assumed 



to be higher because of molecular bands, a higher surface 
temperature at 4.5 \im is interpreted as a sign of a temper- 
ature increasing with height in the atmosphere. Some hot 
Jupiters are thus thought to possess "hot stratospheres", 
likely caused by a visible-light absorber at high altitude. 

Since the brightness temperature observed for 
HD 189733b at 4.5 [im is smaller than at 3.6 (i.m, a 
normal temperature profile (i.e. decreasing with height) is 
inferred. 

However, the "dusty" scenario puts this argument on 
its head. If the opacity difference between the two bands is 
inverted, as suggested by the phase curve results, then the 
relation between the 3.6/4.5 temperature difference and the 
temperature profile is inverted as well. A smaller brightness 
temperature at 4.5 (xm now implies a temperature proflle 
increasing with height, i.e. a "hot stratosphere". In this in- 
terpretation, HD 189733b would have an inverted tempera- 
ture proflle as well. The relation between the temperature 
brightness at 3.6 and 4.5 \im would be caused by an inverted 
ratio of opacities, combined with a rising temperature with 
altitude. 

We note that an inverted temperature proflle is a nat- 
ural consequence of the "dusty" picture. If the opacity de- 
creases from visible to near-infrared, it implies that the in- 
coming flux from the star is absorbed higher than the out- 
going thermal flux is emitted, i.e. an "anti-greenhouse" situ- 
ation. This tends to produce a temperature inversion above 
the infrared photosphere (like the ozone layer on Earth). 
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Figure 13. Spitzer data for the day-side emission spectrum of 
HD 189733b (blue) and best-fit solar-abundance spectrum from 
K12 (black) - adapted from figure 10 in K12. The red line shows 
a blackbody spectrum at 1200 K for the atmosphere of the planet 
(assuming 5000 K for the star). The dotted green line shows the 
expectation from an empirical model based on the transit and 
phase curve data. 



In the transmission spectrum, the opacity in the visi- 
ble is clearly higher than in the infrared, by about an or- 
der of magnitude. If that extends to deeper layers, the anti- 
greenhouse effect will be very strong. 

To first order, the sign of the T-P profile near the pho- 
tosphere for a hot Jupiter is set by the ratio of opacities at 
the wavelengths of the incident starlight to the opacity at 
the wavelengths of the outgoing thermal radiation from the 



planet. In the notation of Hubbard et al. (2001 1, a "Green- 



house factor" 7 can be defined as 7 = Ks/ki where k are 
the fiux-weighted opacities and s and I denote the short- 
wave and longwave opacities. 7 > 1 corresponds to inverted 
temperature profiles. 

This picture is modified, however, if the dust affects 
the opacities by a mixture of scattering and absorption. 



Heng et al. (20121 considered the effect of scattering in the 



simplified analytical framework of Guillot (20101. Then the 



temperature-pressure profile will depend on the absorption 
to scattering ratio. 

We computed the expected T-P profile from the [Heng] 



et al. (20121 relation, assuming that the haze is purely scat- 



tering and that the absorption is dominated by the atomic 
and molecular lines expected from the models. This corre- 
sponds to 7 = 10 and ^ = 0.1 in the notation of | Heng et al.| 
( |2012[ ). ^ parameterises the ratio of scattering to absorption 
in the shortwave. Figure [14] shows the resulting T-P profile. 
A gradient of dT/d(z/H) ~ 200 K is found between the 
longwave and shortwave photospheres (we note that the an- 
alytic approximation predicts an isothermal profile near the 
infrared photosphere, but this is a probably unrealistic fea- 
ture due to the two-band approximation used in separating 
the fiux into a shortwave and longwave component). 

This value is compatible with the measured dayside 
temperature difference between 3.6 and 4.5 \im (136 ± 
14 K) combined with the pressure difference inferred from 
the phase curve amplitude (~ 2), which corresponds to 
dT/d{z/H) ~ 136/(2/e) = 184 K. 

The green line on Fig.[T3]shows the expected secondary 
eclipse spectrum obtained from a model inspired entirely 
from the observations, without assumptions from atmo- 
sphere models. The opacity is assumed to be constant with 
wavelength throughout the infrared, following the transmis- 





=0.0670 cm' q' 
=0.0067 cm' g-' 



Figure 14. Temperature-pressure profile from the Heng et al. 
(2011) analytic model for 70 = 10 and § = 0.10, corresponding to 
a strong anti-greenhouse effect by a highly scattering dust. The 
resulting visible albedo is 0.5, Ks and kq are the shortwave and 
total absorption coefficients. The dotted and dash-dotted lines 
indicate the position of the visible and infrared photospheres. 



sion spectrum, except in the 3.6 \im channel. At 3.6 \im the 
photosphere is assumed to be one scale height higher in pres- 
sure, as suggested by the observed day-night temperature 
contrast. The temperature gradient near the photosphere is 
taken to be -1-200 K per scale height, as given by the analytic 
expression with a scattering haze. We adopted a T=5000 K 
Kurucz spectrum for the star. 

The result is remarkably close to the observed Spitzer 
data, given the lack of tuneable parameters. 

6.2 Pervasive dust 

We therefore converge on the following model for the 
"dusty" scenario of the atmosphere of HD 189733b: 

- Grains are present throughout the transparent part of 
the atmosphere. The mean size of grains decreases towards 
lower pressures, because of the balance between uplifting by 
vertical mixing, and gravitational settling. 

- Scattering by dust dominates the transmission spec- 
trum, causing the transit radius to diminish monotonically 
from 300 nm to 1 \im, and possibly all the way to the mid- 
infrared, masking the pressure-broadened wings of alkali- 
metal absorption and most molecular bands. 

- Dust scattering and absorption in the visible, above 
the level of the thermal-infrared photosphere, causes a tem- 
perature inversion. The corresponding anti-greenhouse effect 
cools the mid atmosphere. 

- The opacity decreases significantly from 3.6 to 4.5 \im 
because of large dust grains, causing the phase curve to be 
more pronounced in the first passband than in the second, 
contrary to expectations if molecules dominate the infrared 
opacities. 

- Because the heat is deposited and transported at lower 
pressure than in the clear-atmosphere case, the radiative 
timescales are lower, and the eastward shift of the hot-spot 
consequently lower. Another observable consequence of the 
same phenomenon is that the day-night amplitude is higher 
than expected from clear-atmosphere models, particularly 
at 3.6 [i.m. 
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Such highly modified opacities compared to a dust-free 
model, especially in the visible, will have significant conse- 
quences for the atmospheric circulation. Overall, the deposi- 
tion of incoming starlight will be moved to lower pressures, 
lowering the radiative timescale. This will tend to make the 
redistribution of heat by atmospheric winds less efficient. 



Perna et al. ( 2012 I have studied the first order effect of vary- 



ing the ratio of visible to infrared opacities on the circulation 
and the eastward drift of the longitude of the hot spot. 

The higher deposition of stellar energy in the atmo- 
sphere will also affect the energy budget of the planet. 
Many hot Jupiters have anomalously large radii, and this 
is thought to be due to a transformation of some of the 
incoming stellar radiation energy into internal entropy in 
the planet. The exact mechanism for this transformation is 
not yet determined. The two leading candidates are (i) deep 
dissipation of the day-night currents (Showman & Guillot 
2002|) (ii) m agnetic drag ( Batygin &: Stevenson|20"T0 Perna 
et al.|20T0 l. If the deposition of incoming starlight is moved 
to lower pressure, the efficiency of both mechanisms could 
drop considerably. 

Obviously, the "dusty" interpretation of HD 189733b 
is at this stage only a possibility. The issue is under- 
constrained by the present data. Given the complexity of 
planetary atmospheres, the patchiness of the observations, 
and the difficulty of controlling instrumental effects in the 
data, a clear-atmosphere interpretation with various addi- 
tional effects remains a valid option. 

What kind of observations could allow us to discrim- 
inate between the two pictures? One prime candidate is 
the observation of the reflection spectrum of the planet in 
the visible and near-infrared. This is challenging but within 
reach of HST. A high albedo in the visible would indicate 
that scattering by condensates affects the zenith geometry 
as well as the grazing incidence. On the contrary, a very 
dark albedo would suggest that the hazes seen in the trans- 
mission spectrum are restricted to a high altitude and do 
not affect the visible opacities near the photosphere. The 
colour dependence of the albedo would provide a clue on 
the importance and grain size distribution of the clouds. 

Given the large amount of observations devoted to this 
planet, and the difficulty of correcting for the variability of 
the parent star, another way forward is to collect similar ob- 
servations for other planets, in the hope that, as happened 
for brown dwarfs, ensemble data will be suggestive of the 
overall effect of condensates if they are a dominant factor 
for some hot Jupiters. Ongoing HST and Spitzer observa- 
tion campaigns will address this issue. A dozen hot Jupiters 
may soon be observed with enough accuracy for a tentative 
examination of the "dusty atmosphere" hypothesis. 



6.3 Origin and composition of the condensates 

Thermal condensates are the most obvious candidate for 
the haze and clouds of HD 189733b. Scattering needs to 
dominate over absorption in the visible for the small, high- 
altitude grain, which favours ^ 0.1 (xm enstatite grains. A 
solar abundance of silicates provides enough grains to ac- 



count for the observed height of the scattering (see Lecave- 
lier Pes Etangs et al.||2008[ ). 

Since the condensation temperature of silicates is above 
1500 K and the photosphere temperatures in HD189733b are 



in the 900-1400 K range, these grains would have to form 
in deeper parts of the planetary envelope, then be trans- 
ported upwards. This is more difficult to conceive in hot 
Jupiters than in brown dwarfs, because the strong stellar 
forcing suppresses convection in the atmosphere. However, 
the day-night recirculation may imply a substantial amount 
of vertical motion. In some models, convection operates in 



the cooler night side (Burkert et al. 2005 Dobbs-Dixon & 



Lin 20081. Gas-phase silicates dredged up from the night 



side could spread to the whole atmosphere, in equilibrium 
with gravitational settling. 

One drawback of this scenario is that silicates are 
strongly absorbent in the mid infrared (beyond 8 vim), but 
the transmission data does not suggest increased absorption 
in the 8 [im and 24 [xm bands. 

Another type of candidate for condensates in hot 
Jupiter atmospheres are photochemical compounds, for in- 
stance sulphur and carbon soots produces by stellar UV ra- 
diation in the upper part of the atmosphere. Because HD 
189733 is a very active star, the planet receives a large 
amount of UV, making photochemistry more likely. Accord- 



ing to Zahnle et al. {2009a b) though, carbon photochemi- 



cal products are expected to be highly absorbent, and sul- 
phur products should produce a large feature shortwards of 
400 nm, neither of which correspond to the observed trans- 
mission spectrum of HD 189733b. 

These two possibilities are not exclusive. Photochemical 
hazes can seed the formation of condensate clouds. This kind 
of interaction is known to operate on Earth and in planetary 
atmospheres in the Solar System. 

6.4 Re-interpreting the two classes of hot Jupiter 
atmospheres 

There is a stark contrast between the two hot Jupiters 
for which extensive spectroscopic data has been ob- 
tained, HD 189733b and HD 209458b. The atmosphere of 
HD 209458 b appears very transparent, with a low albedo, 
and sodium and possibly titanuim oxide absorption in the 
red, and so transparent in the blue that scattering by the Hy- 
drogen molecule becomes detectable. By contrast, the haze 
on HD 189733 b is sufficient to elevate the effective transit 
radius to a few millibars. 

HD 209458b is the prototype of the inverted tempera- 
ture profile hot Jupiters (measured by Ts.g < 74.5), while 
HD 189733b represents the class with Ts.e > T4.5. [Fortneyl 



et al. ( 2008 1 used this dichotomy to separate hot Jupiters 



in two classes, the "pM" and "pL" classes, in analogy with 
the M and L spectral types of brown dwarfs. They spec- 
ulate that in the hotter pM atmosphere, gaseous TiO and 
VO absorbs the visible star light above the infrared photo- 
sphere, thus causing a temperature inversion, whereas in the 
pL classes, the temperature is low enough for TiO and ViO 
to condensate out of the gas phase. 

Subsequent data have not supported a clean correlation 
between atmospheric temperature and the Ta.e/Ti.s ratio. 



Knutson et al. (20101 has noted a strong correlation be- 



tween activity in the host star and apparent temperature 
inversion inferred from the emission spectrum, which could 
correlate the temperature inversion with the UV flux of the 
star instead. 

However, according to our results, it is possible that 
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both HD 189733b and HD 209458b have inverted tempera- 
ture profiles, the difference being that dust on HD 189733b 
modifies the opacities so that the effective temperatures in 



the 3.6 and 4.5 ^.m bands are swapped (see paragraph 6.1 1. 
Could this dichotomy be representative of two important 
classes of hot Jupiters? In that case, the two categories of 
hot Jupiters in Tss/T^.z, ratio would not be caused by a 
temperature inversion, but by an opacity inversion between 
the two bands, relative to clear-atmosphere models. 

HD 209458b and HD 189733b are representative of their 
class in terms of stellar activity as well as temperature. The 
host star HD 209458b is very quiet, while HD 189733 is very 
active. Thus we need spectroscopic data for more object to 
understand if the presence of condensates responds to tem- 
perature or to UV irradiation. A dependence with tempera- 
ture would suggest condensation clouds, whereas a link with 
UV activity would suggest photochemical processes. 

The two categories of hot Jupiters would then corre- 
spond to clear atmospheres around quiet stars, and dusty at- 
mosphere around active stars. This hypothesis can be tested 
on present and future observations of hot gas giant planets 
with Spitzer and HST. 

6.5 Link with brown dwarfs and young Jupiters 

The temperatures in the atmospheres of hot Jupiters is com- 
parable to that of L-T type brown dwarfs. It is thought that 
the colours of L dwarfs are explained by the appearance of 
clouds, which then sink below the photosphere and become 
invisible in T dwarfs. 

Brown dwarf atmospheres have two fundamental differ- 
ences compared to hot Jupiters: the gravity is much higher, 
and the dominant source of energy transfer is convection 
rather than radiation and advection. Nevertheless, the sim- 
ilarities in temperature and composition are sufficient for 
brown dwarf studies to inform our understanding of hot 
Jupiters. 

HD 189733b has a temperature that falls in the T-type 



range (700-1400 K). However, Marley et al. (20121 pointed 



out that, from the point of view of clouds, planets and brown 
dwarf temperatures should not be compared directly. Other 
things being equal, the effect of a smaller gravity allows for 
the persistence of clouds above the photosphere at lower 
temperature. The fundamental reason is that, at lower grav- 
ity, a given pressure corresponds to more mass. Because the 
opacity of grains is independent of gas pressure, while the 
opacity of atomic and molecular lines grows with pressure, 
cloud opacity will remain important longer at lower gravity. 

Figure [15] shows the position of some brown dwarfs, hot 
Jupiters and directly imaged planets in a temprature-gravity 
diagram. The dashed line shows the gravity dependence of 
the disappearance of clouds below the photosphere accord- 



ing to Marley et al. ( 2012 1. Grey zones and grey symbols in- 
dicate the possibly "dusty" objects, i.e. atmospheres where 
haze/cloud opacities seem required to explain the observed 
spectra. The atmosphere of the wide-orbit, young planets 



found by direct imaging around HD 8799 (e.g. Marois et al. 
|2008[ ) have temperatures similar to HD 189733b. The pres- 
ence of clouds has been inferred in their atmosphere from 
their infrared colours ( Currie et al.|20lT I . 

This provides another angle to tackle the behaviour of 
clouds in the gravity-irradiation plane of sub-stellar objects. 



In spite of the fundamental difference in irradiation regime 
and gravity, the observations suggest the default hypothesis 
that the same kind of clouds form in hot Jupiters than in 
less irradiated young planets of the same temperature, and 
hotter brown dwarfs. 

As more planets and brown dwarfs are characterised 
by direct-imaging observations, and more spectroscopic data 
are collected on hot Jupiters and cooler transiting Jupiters, 
we may hope to obtain a fuller picture of this hot substel- 
lar class of atmospheres and the role of high-temperature 
clouds. 



T SUMMARY AND CONCLUSION 

In this paper we have attempted to build a complete trans- 
mission spectrum of the atmosphere of the hot Jupiter 
HD 189733b from the UV to the infrared. We have used 
three different instrument on HST to cover the range from 
300 nm to 1 |i.m. Other HST and Spttzer observations extend 
the coverage to 24 (j.m in a sparser manner. 

The presence of star spots on the surface of the host star 
HD 189733, and of instrumental systematics in the space 
data at the high level of precision required, dominate the 
error budget. We have devoted much care to accounting for 
these two factors. In both cases, we have used Gaussian pro- 
cesses to model the nuisance factors in a non-parametric and 
Bayesian way. These approach yields larger but, we believe, 
more realistic uncertainties than the usual parametric fits 
to both the instrumental systematics and the star spot vari- 
ability. 

We find that the atmospheric transmission spectrum of 
HD 189733b is globally featureless from 300 nm to 1 \im, ris- 
ing towards the blue with a slope compatible with Rayleigh 
scattering by small (< O.ljxm) grains of condensates. The 
strong cores of the sodium and (possibly) potassium dou- 
blets are the only features rising above the continuum. 

In the infrared, 1-24 \im, we find that the uncertain- 
ties preclude any definite conclusion about the shape of the 
transmission spectrum. The data is compatible with an ex- 
tension of the featureless spectrum in the visible, with a 
fiattening of the slope, or with weak molecular features. It 
is not compatible with the default, clear-atmosphere, solar- 
abundance models. 

Overall the transmission spectrum suggests an extended 
presence of haze/clouds in the atmosphere of the planet. 

We then combine the information from the transmission 
spectrum with the indications from the day-side emission 
spectrum and the phase curves measured with Spitzer in the 
infrared. We find that several of the anomalies in these data 
could be explained by the prevalence of condensates in the 
atmosphere. Notably, the high amplitude of the phase curve 
at 3.6 nm is difficult to reconcile with clear-atmosphere mod- 
els, and emerges more naturally in the presence of conden- 
sates. The lower-than-expected eastward shift of the hottest 
point on the day side would also be a natural consequence 
of the presence of clouds. 

Overall, while several interpretations remain possible, 
we find that the present data suggest the possibility that 
opacity from condensates dominates the atmosphere of HD 
189733b, with important consequence not only on the trans- 
mission and emission spectrum, but also on the atmospheric 
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Figure 15. Brown dwarfs and planets in a temperature-gravity diagram, showing tlie approximate position in temperature of the four 
classes of brown dwarf spectra, M, L, T and Y. Grey areas show the domains where the infrared colours are believed to indicate the 
presence of clouds, i.e. L type and possibly Y type. Ellipses show the position of the direct-imaging planets HD 8799c, d, HD 189733b, 
and two other well-studied hot Jupiters, HD209458b and WASP-12b. Grey ellipses indicate objects possibly dominated by clouds. The 
gravity dependence of cloud effects on brown dwarf and planet s is from|Marley et al.| | |2012[ l. The dotted line shows the M = 13Mj 
planet/brown dwarf limit for R = IRj. Data for HR 879 9 from \Marley et al.\ \201^ , for HD 189733b from \Knutson et al.\ \201Sl^ , for 
HD 209458b and WASP- 12b from \Cowan & Agol{ 1(2011!^ . 



structure, circulation and evolution. In particular, the di- 
chotomy between hot Jupiters with "hot stratospheres" (i.e. 
temperature inversion near the photosphere) and without, 
may in that case be instead a dichotomy of hot Jupiters with 
and without dust, both classes having a hot stratosphere. 

Placing HD 189733b in the wider context of gas gi- 
ant planets and brown dwarfs, we suggest two possibilities: 
(1) the two classes of hot Jupiters corresponds to clear- 
atmosphere planets and to planet in which a high UV flux 
from the active parent star trigger the formation of a photo- 
chemical haze and clouds (2) hot Jupiters like HD 189733b 
have the same silicate/iron clouds as L-type brown dwarfs, 
and as young planets found by direct imaging like the com- 
panions of HR 8799. 

Both hypothesis are testable in the near future by new 
observations. It is necessary to stress that the interpreta- 
tion of the data on HD 189733b presented here is tentative, 
and is not unique. Observational uncertainties are large and 
sometimes poorly constrained, some observations are only 
marginally compatible with others, and the models have 
many parameters since they have to account for global fea- 



tures of the whole planet. A large amount of space-telescope 
time has been devoted to HD 189733b, and it is probable 
that definite progress will have to wait for a new generation 
of observatories, such as JWST or dedicated space projects. 
In the meantime, more insight may come from the character- 
isation of other comparable planets. A few more hot Jupiters 
are accessible to spectroscopic characterisation with HST, 
and new young planets are expected to be discovered in the 
coming years by direct-imaging programmes. 

Our results also have implications for the study of exo- 
planet spectra, and the search for biomarkers on terrestrial 
exoplanets in the longer terms. By necessity, information is 
often inferred for exoplanets by fitting suites of model spec- 
tra in a Bayesian fashion to a few Spitzer passbands, some- 
times complemented by a few ground-based measurements 
at shorter wavelength. These "spectral retrieval" suites gen- 
erally do not include condensates and clouds (which would 
make the Bayesian integration intractable and too depen- 
dent on prior assumptions about the properties of the 
clouds). If our experience with HD 189733b is any guide, 
we caution against taking these results too seriously. The 
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results of Lee et al. (20121 for instance are invalidated by 



the contribution of haze or clouds. The implication of the 
case of HD 189733b could be summed up as: beware of in- 
complete spectra. It is clear that fitting a suite of synthetic 
spectra to a restricted subset of the data considered in this 
paper would lead to very misleading conclusions. Each time 
more extended data has been forthcoming, it has flatly de- 
feated the expectation. 

This has implications for the design of instruments and 
space missions for the study of exoplanet atmospheres. It 
would tilt the balance towards an extended spectral cov- 
erage, not too narrowly focussed on expected features and 
model predictions, and keeping the possible presence of con- 
densates in mind when calculating detection capabilities. 
This may also extend to the search for biomarkers in habit- 
able exoplanets. 

Our results suggest several lines of investigation for hot 
Jupiter atmosphere models. 3-D circulation models can trace 
the formation and circulation of dust in a hot Jupiter like 
HD 189733b, and test whether the distribution implied by 
the transmission spectrum can result from silicate grains. 
Atmosphere structure models can add haze/clouds opaci- 
ties to the radiative transfer terms and study the conse- 
quences on the observable signatures, such as phase curve, 
photosphere temperature and wind speed. Models of pho- 
tochemical pathways can explore the possibility that small, 
transparent grains can be produced by stellar UVs in the 
upper atmosphere of hot Jupiters. 
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APPENDIX A: GAUSSIAN PROCESS 
INTERPOLATION OF THE STELLAR FLUX 

In this Appendix, we explain how we used ground-based 
observations of the brightness of HD 189773, coUected over 
several years, to evaluate the impact of unocculted star-spots 
on transmission spectra obtained during the same period 
with HST and Spitzer. Our goal is to estimate the fraction 
of the stellar disk that is covered in spots at the time of each 
transit observation, and the resulting dimming of the star 
(relative to the unspotted flux, which is not known) in the 
relevant bandpass. 

We proceed as follows: 

(i) combine long-term out-of-transit monitoring of 
HD 189733 from several observatories in different band- 
passes; 

(ii) model the resulting time-series using a quasi-periodic 
Gaussian process (GP) model; 

(iii) estimate the star's brightness at the time of each 
transit observation, relative to a putative un-spotted bright- 
ness; 

(iv) convert this to the bandpass of each observation; 

(v) use the results to convert the measured planet-to-star 
radius ratio in each band. 

There are 2 free parameters in this process: the spot tem- 
perature, used in steps (i) and (iv), and the un-spotted flux, 
used in step (iii). 



Combining datasets 

Our main source of out-of-transit data for HD 189733 is the 
long-term monitoring program carried out with the 0.8-m 
Automated Patrol Telescope (APT) at Fairborn Observa- 
tory ( [Henry fc Winn||2008[ ). The APT uses photomultiplier 
tubes to gather photometry in the Stromgren b and y filters 
simultaneously, resulting in an effective 'b+y' bandpass. This 
program started in 2005 and is ongoing, with two-month- 
long breaks each semester (in August and January). There 
are typically one or two points per night, resulting in a total 
of 868 observations. 

We also used observations taken with the 40-inch tele- 
scope at Wise Observatory in November 2009 using the 



Bessel R filter (Sing et al. 20111. There were 24 observa- 



tions, with a typical sampling of 1 point per night. To ac- 
count for the difference in bandpass with the APT data, the 
Wise data were scaled by a factor 1.12. This factor is simply 
Rfj^'^iS, T), as defined in Section 2.2.4, where S = 0.01, i and 
j corresponding to the Wise R and APT b + y bandpasses 
respectively, and we adopted T = 4250 K for the spots. Fur- 
thermore, a constant offset was added to the Wise data, to 
bring its median magnitude level with that of the APT data 
in the two nearest semesters. 



GP regression basics 

We model the out-of-transit brightness variations of 
HD 189733 as a GP, which enables us to interpolate these 
measurements to the times of the HST and Spitzer observa- 
tions, subject to certain smoothness constraints. Specifically, 
we use a quasi-periodic kernel to account for the fact that 
the spot-induced variability has a dominant periodicity (the 



rotation period of the star), but its amplitude and phase 
evolve. An important feature of the GP model is that it pre- 
dicts a probability distribution for the star's magnitude at 
any specific date, i.e. our interpolates are associated with 
robust error bars. 

Here we summarise the main steps in GP regression 
very succinctly. The interested reader will find an excel- 
lent textbook-level introduction to GPs in IRasmussen fcl 



Williams { 2006 1 , and a more detailed description of the GP 



regression process as applied to astrophysical time-series in 



the appendix of Gibson et al. (2012b I 



Under a GP model, the joint probability distribution of 
any finite set of N outputs (in the present example, out-of- 
transit brightness measurements for HD 189733) is assumed 
to be a multi-variate Gaussian distribution. This distribu- 
tion is fully specified by an A'^-element vector of mean values, 
and an N X N covariance matrix. In the present application, 
the mean vector was assumed to be constanl[f]and the covari- 
ance matrix is used to describe both the intrinsic brightness 
variations of the star and the observational noise. 

The elements of the covariance matrix are specified by 
a covariance function or kernel, which takes a pair of inputs 
(in the present case, the observation times) and returns the 
covariance between the corresponding pair of outputs. This 
function is controlled by a number of parameters, which are 
known as the hyper-parameters of the GP, since they spec- 
ify the covariance properties of the GP, rather than directly 
specifying an output value at a particular time. Together, 
the kernel function and hyper-parameters specify a prior 
distribution over random functions sharing the same covari- 
ance properties. Once this prior is conditioned on the data, 
it gives a probability distribution for the expected output (s) 
at any given set of input (s), which is also Gaussian. The 
choice of kernel function and of hyper-parameters are clearly 
critical, and need to refiect whatever information we have 
about the underlying physical process. As we will see be- 
low, there are many different ways of comparing kernels and 
hyper-parameters, ranging from a fully Bayesian to maxi- 
mum likelihood and other, more ad-hoc methods. 



Constructing a kernel for star spot variability 

We seek a kernel which appropriately describe the kind of 
brightness variations arising from a time-evolving distribu- 
tion of spots on a rotating stellar surface. We take the fol- 
lowing commonly used kernels as our basic building blocks. 
White noise with standard deviation a is represented by: 



kwN{t, t') — cr^I, 



(Al) 



where I is the identity matrix. The squared exponential (SE) 
kernel is given by: 

(t^t'r 



ksE{t,t') = A exp 



2P 



(A2) 



where A is an amplitude and I is a length scale. This gives 
rather smooth variations with a typical time-scale of I and 
r.m.s. amplitude A. The rational quadratic (RQ) kernel is 
given by: 

^ In fact, the time-series was median-subtracted before mod- 
elling, and the mean vector was set to zero. 
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kR,Q{t,t') 



jt-t'f 

2aP 



(A3) 



where a is known as the index. IRasmussen fc WiUiams] 
(2006) show that this is equivalent to a scale mixture of 
SE kernels with different length scales, distributed accord- 
ing to a Beta distribution with parameters a and This 
gives variations with a range of time-scales, the distribution 
peaking around I but extending to significantly longer peri- 
ods (but remaining rather smooth). When a — > oo, the RQ 
reduces to the SE with length scale /. There are many more 
types of covariance functions in use, including some (such 
as the Matern family) which are better suited to model 
rougher, less smooth variations. However, the SE and RQ 
kernels already offer a great degree of freedom with rela- 
tively few hyper-parameters, and covariance functions based 
on these are sufficient to model the data of interest satisfac- 
torily. 

A periodic covariance function can be constructed from 
any kernel involving the squared distance [t — t')^ by replac- 
ing the latter with sin^[7r(t — t')/P], where P is the period. 
For example, the following: 



kpcT,SE{t,t') = A exp 



Ht-t')/p] 

2L2 



(A4) 



gives periodic variations which closely resemble samples 
drawn from a squared exponential GP within a given kernel. 
The length scale L is now relative to the period, and letting 
L — > oo gives sinusoidal variations, whilst increasingly small 
values of L give periodic variations with increasingly com- 
plex harmonic content. Similar periodic functions could be 
constructed from any kernel. Other periodic functions could 
also be used, so long as they give rise to a symmetric, posi- 
tive definite covariance matrix - sin'^ is merely the simplest. 
As described by Rasmussen & Williams' (2006), valid 



covariance functions can be constructed by adding or multi- 
plying simpler covariance functions. Thus, we can obtain a 
quasi-periodic kernel simply by multiplying a periodic ker- 
nel with one of the basic kernels described above. The latter 
then specifies the rate of evolution of periodic signal. For 
example, we can multiply equation | A4| with yet another SE 
kernel: 



kQP,SE{t,t') 



A exp 



sm^[TT{t-t')/P] {t-t'f 



2L2 



2P 



(A5) 



to model a quasi-periodic signal with a single evolutionary 
time-scale I. In the case of activity-induced stellar brightness 
variations, which are caused by the evolution and rotational 
modulation of active regions, one may expect a range of both 
periodic covariance scales L and evolutionary time-scales /, 
corresponding to different active region sizes and life-times 
respectively. This can be achieved by replacing one or both 
of the SE kernels in equation |A5| by RQ kernels. 

Finally, we can also allow for correlated noise on short 
to moderate time-scales by including a separate, additive SE 
or RQ kernel. 



Training the GP and comparing kernels 

The marginal likelihood, which is the product of the predic- 
tive probabilities for the observed outputs, provides a 'good- 
ness of fit' measure which can be used to optimise the hyper- 



parameters for a given kernel (this is known as training the 
GP). When doing this, we exploited the fact that simple an- 
alytical expressions exist for the derivatives of the marginal 
likelihood with respect to the hyper-parameters to speed up 
the optimisation using conjugate gradient methods. 

However, the marginal likelihood surface can be rather 
complex, specially for (quasi-)periodic kernels. This means 
that there is a tendency for the optimiser used to be- 
come trapped in local optima. When studying a partic- 
ular dataset, as in the present case, this is circumvented 
by (a) carefully choosing the initial guesses for the hyper- 
parameters, based on a visual inspection of the data, and 
(b) repeating the GP training process using different initial 
guesses. 

Ideally, we would prefer to map out the entire marginal 
likelihood surface, to ensure that we have found the global 
optimum, and to enable us to marginalise with respect 
to the hyper-parameters in order to compare different 
kernels. There are several ways of doing this: evaluating 
the marginal likelihood at a grid of points in the hyper- 
parameter space, using global optimisers, Markov Chain 
Monte Carlo (MCMC) methods, or Bayesian quadrature 
(which consists in modelling the likelihood surface itself, 
for example using a GP, in order to evaluate the integrals 
involved in the marginalisation process). Each evaluation 
of the marginal likelihood is computationally costly, mak- 
ing MCMC approaches prohibitive. We experimented with 
global optimisers, but with limited success. The Bayesian 
quadrature option is certainly the most promising, but a 
significant amount of code development would be needed to 
implement the necessary nested hierarchical models. There- 
fore, the results presented in this document are based on grid 
sampling around a manually selected initial guess, combined 
with a local optimiser. 

We used two methods to compare different kernels: 
comparing the marginal likelihood obtained with the best- 
fit hyper-parameters in each case, and leave-one-out cross- 
validation (LOO-CV). Cross-validation consists in training 
the model on a subset of the data and then testing its ability 
to predict the remaining subset. In the case of GP, this is 
done by measuring the likelihood of the test set given the 
partially trained GP, as given by the predictive distribution 
for that subset. This process is carried out repeatedly, ex- 
cluding different subsets in turn, and multiplying the results 
together to obtain a 'pseudo-likelihood'. Rather than mea- 
suring the 'goodness of fit', as the marginal likelihood does, 
this pseudo-likelihood measures the predictive ability of the 
model. Given that our task is making predictions at times 
where we do not have observations, this seems an appro- 
priate way of comparing kernels. LOO-CV is a special case 
where each data point is 'left-out' in turn. This is generally 
prohibitive computationally, as the model must be trained 
anew for each data point that is excluded. However, in the 
case of CPs there is a neat shortcut that allows the pseudo- 
likelihood for LOO-CV to be computed directly from the 
inverse of the full covariance matrix, making it a workable 
proposition. 

Choosing a model for HD 189733 

We modelled the HD 189733 data with all the individual 
kernels detailed above, as well as a number of different com- 
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binations. We used LOO-CV to compare different kernels, 
but we also performed a careful visual examination of the 
results, generating predictive distributions over the entire 
monitoring period and over individual seasons. Interestingly, 
we note that, in practice, LOO-CV systematically favours 
the simplest kernel which appears to give good results by 
eye. 

We experimented with various combinations of square- 
exponential (SE) and rational-quadratic (RQ) kernels to 
form quasi-periodic models, and found that using an RQ 
kernel for the evolutionary term significantly improves the 
LOO-CV pseudo-likelihood relative to an SE-based evolu- 
tionary term. It makes very little difference to the mean of 
the predictive distribution when the observations are well- 
sampled, but it seems to increase the predictive power of the 
GP away from observations, as it allows for a small amount 
of covariance on long time-scales, even if the dominant evo- 
lution time-scale is relatively short. On the other hand, we 
were not able to distinguish between SE and RQ kernels for 
the periodic term (the two give virtually identical best-fit 
marginal likelihoods and pseudo-likelihoods), and therefore 
opted for the simpler SE kernel. 

To describe the observational noise, we experimented 
with a separate, additive SE kernel as well as a white noise 
term. However, we found that this did not significantly im- 
prove the marginal or pseudo-likelihood, and the best-fit 
length-scale was comparable to the typical interval between 
consecutive data points. We therefore reverted to a white- 
noise term only. The final kernel was thus: 



od(t,t')= A^expl- 



!in^[Tr(i-t')/P] 



1 4- 



(A6) 



The best-fit hyperparameters were A — 6.68 mmag, 
P = 11.86 days, L = 0.91, a = 0.23, I = 17.80 days, and 
a — 2.1 mmag. Our period estimate is in excellent agree- 



ment with Henry & Winn ( 2008 1. We note that very similar 



best-fit hyper-parameters were obtained with the other ker- 
nels we tried (where those kernel shared equivalent hyper- 
parameters). The relatively long periodic length-scale L in- 
dicates that the variations are dominated by fairly large ac- 
tive regions. The evolutionary term has a relatively short 
time-scale I (about 1.5 times the period) but a relatively 
shallow index a, which is consistent with the notion that in- 
dividual active regions evolve relatively fast, but that there 
are preferentially active longitudes where active regions tend 
to be located (as inferred from better-sampled long-duration 
CoRoT light curves for similarly active stars). 



this predicted difference. This can be obtained directly from 
the covariance matrix of the predictive distribution: 

"'v2\vi = cov(xi, xi) + cov(x2, 2:2) — 2cov(a;i, a;2). (A7) 

The results are reported in Table [S] We take as our ref- 
erence the maximum fiux predicted by the GP at any time 
throughout the 5-1- years of monitoring with the APT. As- 
suming that the spots are dark, this represents a lower limit 
on the spot-free flux. The flux differences relative to this level 
can readily be converted to differences relative to another, 
arbitrarily higher fiux by adding a constant to all of them. 
Note that we worked in magnitude units when training and 
conditioning the GP, but have been careful to convert the 
results to fiux units before tabulating them. We computed 
the fiux drop values assuming a spot temperature of 4000 K. 
Changing the spot temperature changes the training data 
very slightly, because it changes the conversion factor be- 
tween the Wise and APT bandpasses. However, this makes 
only a very small difference to the fiux estimates for transit 
observations occurring during the same season as the Wise 
observations, and none at all for the other transit observa- 
tions, so we have not tabulated the estimated fiux differences 
assuming other spot temperatures. This does not preclude 
trying out various spot temperatures when converting the 
fiux drops reported in Table [3] to the bandpasses of the var- 
ious transit observations. 



APPENDIX B: UNCERTAINTIES ON PHASE 
CURVE AMPLITUDES DUE TO STELLAR 
ACTIVITY 



In [Knutson et al.| ( |2012[ ) , the quasi-periodic GP trained on 
the ground-based APT data (as described in Appendix A) 
was used to predict whether the stellar flux during the 
Spitzer observations was linear or underwent an inflection, 
and if so what direction the inflection went. An MCMC code 
was then used to adjust the coefficients of a model that con- 
tained these terms for the stellar variability, plus a series of 
systematics terms, plus planetary phase-curve terms. The 
final uncertainties on the phase curve amplitudes were de- 
rived directly from these MCMC runs. 
Here we try to address two questions: 

(i) Is the predictive capability of the GP model sufficient 
to use it in that way? 

(ii) To what extent can HD 189733's variability mimic a 
phase curve? How does that compare with the errors on the 
amplitude reported in the K12 paper? 



Results 

Once the GP has been trained and conditioned on the avail- 
able data, we can compute a predictive distribution for any 
desired set of input times. This predictive distribution is a 
multi-variate Gaussian, specified by a mean vector and a 
covariance matrix (see e.g. [Gibson et al.|2012b[ for details). 
Here we are interested in estimating the difference between 
the predicted fiux at two input times, which is just the differ- 
ence between the corresponding elements of the mean vector. 
We also want an estimate of the uncertainty associated with 



Predictive power of the GP 

As shown in Fig. |B1| reproduced from K12, the Spitzer ob- 
servations at 3.6 and 4.5 [im each occur after the end of 
a ground-based monitoring season. As a result, the CP's 
ability to predict the overall level, amplitude and phase of 
the quasi-periodic modulation of the stellar fiux during the 
Spitzer observations is somewhat limited. This is illustrated 
by Fig. |B2[ where we have drawn random samples from the 
GP around the time of the Spitzer observations. 

Nonetheless, this exercise can be used to check what 
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Figure Bl. Timing of the Spitzer pliase curve observations at 
3.6 (left) and 4.5 ^.m (right), indicated by the vertical blue hues, 
compared to the APT observations (black points). The mean and 
1-(T confidence interval of the predictive distribution of the GP 
model as shown by the red line and pink shaded area, respectively. 
Reproduced from K12. 




Figure B3. MOST 2007 observations of HD 189733 (green dots) 
together with the predictive mean (red line) and standard devia- 
tion (pink shaded area) of the GP model trained on the APT data 
alone, but conditioned on the MOST data. Note that the stan- 
dard deviation includes the white noise level of the APT data, 
which is significantly larger than that of the MOST data. 




Figure B2. Random realisations of the GP (coloured lines) 
drawn around the time of the Spitzer observations. 



fraction of the time the GP samples are better fit by a lin- 
ear or quadratic function, and to check if the quadratic term 
tends to be positive or negative. This is what was done in 
K12 in the case of the 3.6-vim passband, the GP samples 
tend to undergo a minimum, so K12 used a quadratic term 
and forced it to be positive. In the case of the 4.5-nm pass- 
band, the GP samples tend to decrease monotonically, so 
K12 used only a linear term. 

One Umitation of the approach above is that the sam- 
ples are still drawn from a single GP, whose covariance 
hyper-parameters were trained on the APT data. The un- 
certainty on these hyper-parameters is not accounted for. In 
particular, the rather sparse sampling of the APT data (one 
point every day or two during the observing seasons) means 
that we cannot constrain the periodic length scale well. The 
best-fit GP has a fairly large periodic length scale, making it 
rather smooth, with only one or two inflections per period. 
But if there were more short-term variations in a period, 
we would not be able to distinguish them from white noise 
using the APT data alone. Since the phase curve signal op- 
erates on the time scale of the planet orbital period (2.2 
days), these short-term variations would affect the phase 
curve measurements. 

We can use the 2007 observations of HD 189733 by 
the MOST space mission (kindly provided by Brice KroU), 
which have much better precision and time sampling than 
the APT data, and were taken outside the APT observ- 
ing season, to test the importance of the limitations noted 
above. The MOST data provides continuous coverage of the 
brightness of HD 189733 over one month. The transits of 
HD 189733b had already been removed from these data, and 
a linear trend of order 10"* mag per day, suspected to be 
of instrumental origin, had been subtracted. We scaled the 



MOST magnitude variations by a factor 1/0.987 to account 
for the difference between the MOST and APT passbands 
(the former was crudely approximated by a top-hat through- 
put function from 350 to 700 nm), and added a constant to 
the MOST dataset to bring its median magnitude level with 
that of the surrounding two semesters of APT observations. 
Finally, we binned up the MOST data (from ~1 point per 
hour to ~2 points per day, resulting in a total of 69 data 
points). 



In Fig. B3 we have trained the GP (adjusted its hyper- 
parameters) on the APT data, but then conditioned it (ad- 
justed the mean and covariance, but keeping the hyper- 
parameters fixed) on the MOST data (shown in green, after 
removing the transits, and binning to two points per day for 
clarity). The GP is clearly unable to follow the more rapid 
variations in the MOST flux. This shows that the GP trained 
on the APT data does not capture the short-timescale vari- 
ability. In fact the variability it fails to capture is precisely 
on the 1-2 day timescales relevant to the phase curve mea- 
surements. Additionally, we also performed a more stringent 
test, where the GP was trained and conditioned only on the 
APT data, and compared to the MOST observations with- 
out conditioning it on the latter. The results are shown in 
Fig. |B4| They show very clearly that, although the gap be- 
tween the APT and MOST observations is only around one 
stellar rotation period (much like the gap between the APT 
and Spitzer observations) , this is enough for the GP to lose 
its ability to predict the phase and, to a lesser extent, the 
amplitude of the signal. In the light of this new test, we are 
compelled to question the validity of using the GP to con- 
strain the Spitzer flts, even in the very conservative fashion 
adopted in K12. 



Estimating the activity-induced error on the phase 
curve amplitude: 

We know that stellar variability contributes to the uncer- 
tainty on the phase-curve amplitude, but we have shown, 
using the MOST data, that the GP cannot be used to correct 
for this effect reliably. On the other hand, the same MOST 
observations can be used to evaluate the contribution of ac- 
tivity to the uncertainty on the phase curve amplitude, again 
making use of the MOST data. Stellar variability is expected 
to be by far the dominant contribution to the brightness 
fluctuations observed with MOST; any orbital phase varia- 
tions of planetary origin (arising, from example, from star 
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Figure B4. MOST (green dots) and APT (black dots) obser- 
vations of HD 189733 (green dots) together witli tlie predictive 
mean (red line) and standard deviation (pink shaded area) of the 
GP model trained and conditioned on the APT data alone. 



to the Spitzer channel 1 and 2 bandpasses respectively. To 
summarise, stellar variability induces an uncertainty on the 
Spitzer phase curve amplitudes of at least 0.027 and 0.025% 
at 3.6 and 4.5 i^m, respectively. For comparison, the phase 
curve amplitude and uncertainties adopted in K12, derived 
from constrained MCMC fits to the Spitzer light curves, were 
0.124 ± 0.006 and 0.098 ± 0.009% respectively. 

To understand the implications for this in terms of at- 
mospheric physics, we need to convert the phase curve am- 
plitudes to temperature contrasts. This gives uncertainties of 
150 and 78 K for the brightness temperatures at 3.6 and 4.5 
[im, respectively. Thus the new value for the day-side bright- 
ness temperatures measured during the secondary eclipses 
become: 1328±150K at 3.6 \im, and 1192±78K at 4.5 ]xm. 
Our revised error estimates do not bring the detection of 
a phase modulation into question (both remain > 4cr), but 
it does substantially decrease the significance level of the 
difference between the day-night contrast in the two pass- 
bands. 
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Figure B5. Histograms of the amplitude of the phase-curve com- 
ponent of the fitted model (in parts per thousand). The two 
colours correspond to linear and quadratic models for the vari- 
ability terms. The x-axis is in parts per thousand. 



light refiected by the planet) are expected to be tiny, and 
certainly much smaller than the infrared phase-curve ampli- 
tudes. Thus, if all we observe is stellar variability, but we 
try to model is as a combination of a planetary and a stel- 
lar signal, as done in K12, what phase curve amplitudes do 
we obtain? This can be used to estimate the contribution 
from stellar variability to the error on the phase curve am- 
plitude (after applying a suitable scaling for the difference 
in bandpasses). 

To perform this test, we took the (un-binned) MOST 
time-series, extracted sections lasting 2.7 days (approximate 
duration of the Spitzer phase curve observations) and fit 
them with a model containing a polynomial (order 1 or 2) 
to represent the stellar variability, and a sum of sines and 
cosines at the period and half the period of the planet, to 
represent a potential planetary phase curve. The phase of 
the sinusoidal terms was always fixed so that the transit 
would have occurred in the middle of the data segment used 
(as in the Spiteer observations). Figure B5 shows histograms 



of the amplitude of the phase-curve component of the fitted 
model. Regardless of the order of the polynomial term, the 
stellar variability can mimic phase curve signals at the level 
of at least 0.1% in the optical. This needs to be divided 
by a factor of 3.7 and 3.9 to convert it from the MOST 



APPENDIX C: SLOPE OF THE 
TRANSMISSION SPECTRUM WITH 
SETTLING GRAINS 

We calculate the slope of the transmission spectrum, in other 
words the dependence of the effective transit radius r as a 
function of wavelength A, when the opacity is dominated by 
Rayleigh scattering from grains of condensates of effective 
size a, with the density and maximum size of the grains vary- 



ing with pressure ("settling grains" case in Section 4.2.21. 

We use the usual approximation that the effective tran- 
sit radius corresponds to the height in the atmosphere at 
which the opacity reaches a critical value Tcrit. 

If ^(2) is the abundance of grains at height z in the 
atmosphere, then 



r(z, A) oc / n{z)a{a,X)da 



p{z)£^{z)a{a, X)da 



where n is the number density of grains, p the atmospheric 
density and a the cross-section of individual grains. 
For Rayleigh scattering, 

a(a, A) oc a^A~* 

Because of the very steep dependence of the cross- 
section on grain size, the largest grains at a given height 
dominate. We therefore use: 

a{z,X) oc amax(2)A'' 

thus 

Tcrit OC p(z)^{z) {z)\-' 

In the ideal-gas, isothermal approximation, 
p{z) cc p{z) 
and 

p(z) OC e ' 

where p is the gas pressure and H is the atmospheric scale 
height. 

For the scenarios in Section [2] we assume a,nax(2) oc 
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p{z) (settling grains in the molecular regime), and use the 
following dependence of amax on p: 

) oc 

with j3 — corresponding to a flat distribution of grain 
sizes, and /3 = —3 corresponding to an equal partition of 
mass across all grain sizes. This gives 

= p{z)p(zfp{zfX-^=p^'+^^\-'' 

^ g-(7+/3)z/H;^-4 

Taking the logarithm gives: 
4 

z = i^lnA 

7 + 13 

Therefore the slope of the transmission spectrum, in 
units of atmospheric scale height H per In A, is —4/7 in 
the case of flat grain size distribution (/3 — 0), and —1 in 
the case of constant mass fraction distribution of grain sizes 
(/? = -3). 

This paper has been typeset from a TgX/ LM^rjX file prepared 
by the author. 
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